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.18.0 lcov report (development 29732-95c6201d93) Lines: 2065 2174 95.0 %
Date: 2024-11-21 09:08:54 Functions: 227 238 95.4 %
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    40230586 : get_tab(GEN nf, long *N)
      33             : {
      34    40230586 :   GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
      35    40230586 :   *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  1087478756 : _mulii(GEN x, GEN y) {
      41  1757159163 :   return is_pm1(x)? (signe(x) < 0)? negi(y): y
      42  1757017268 :                   : mulii(x, y);
      43             : }
      44             : 
      45             : GEN
      46       22203 : tablemul_ei_ej(GEN M, long i, long j)
      47             : {
      48             :   long N;
      49       22203 :   GEN tab = get_tab(M, &N);
      50       22203 :   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       11557 : tablemul_ei(GEN M, GEN x, long i)
      58             : {
      59             :   long j, k, N;
      60             :   GEN v, tab;
      61             : 
      62       11557 :   if (i==1) return gcopy(x);
      63       11557 :   tab = get_tab(M, &N);
      64       11557 :   if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
      65       11557 :   tab += (i-1)*N; v = cgetg(N+1,t_COL);
      66             :   /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
      67       78491 :   for (k=1; k<=N; k++)
      68             :   {
      69       66934 :     pari_sp av = avma;
      70       66934 :     GEN s = gen_0;
      71      473214 :     for (j=1; j<=N; j++)
      72             :     {
      73      406280 :       GEN c = gcoeff(tab,k,j);
      74      406280 :       if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
      75             :     }
      76       66934 :     gel(v,k) = gerepileupto(av,s);
      77             :   }
      78       11557 :   return v;
      79             : }
      80             : /* as tablemul_ei, assume x a ZV of correct length */
      81             : GEN
      82    23969067 : zk_ei_mul(GEN nf, GEN x, long i)
      83             : {
      84             :   long j, k, N;
      85             :   GEN v, tab;
      86             : 
      87    23969067 :   if (i==1) return ZC_copy(x);
      88    23969067 :   tab = get_tab(nf, &N); tab += (i-1)*N;
      89    23967598 :   v = cgetg(N+1,t_COL);
      90   169916894 :   for (k=1; k<=N; k++)
      91             :   {
      92   145951708 :     pari_sp av = avma;
      93   145951708 :     GEN s = gen_0;
      94  2143476409 :     for (j=1; j<=N; j++)
      95             :     {
      96  1997681548 :       GEN c = gcoeff(tab,k,j);
      97  1997681548 :       if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
      98             :     }
      99   145794861 :     gel(v,k) = gerepileuptoint(av, s);
     100             :   }
     101    23965186 :   return v;
     102             : }
     103             : 
     104             : /* table of multiplication by wi in R[w1,..., wN] */
     105             : GEN
     106       39293 : ei_multable(GEN TAB, long i)
     107             : {
     108             :   long k,N;
     109       39293 :   GEN m, tab = get_tab(TAB, &N);
     110       39293 :   tab += (i-1)*N;
     111       39293 :   m = cgetg(N+1,t_MAT);
     112      154061 :   for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
     113       39293 :   return m;
     114             : }
     115             : 
     116             : GEN
     117    10843934 : zk_multable(GEN nf, GEN x)
     118             : {
     119    10843934 :   long i, l = lg(x);
     120    10843934 :   GEN mul = cgetg(l,t_MAT);
     121    10843910 :   gel(mul,1) = x; /* assume w_1 = 1 */
     122    34447892 :   for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
     123    10840248 :   return mul;
     124             : }
     125             : GEN
     126        2751 : multable(GEN M, GEN x)
     127             : {
     128             :   long i, N;
     129             :   GEN mul;
     130        2751 :   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     5001615 : zk_scalar_or_multable(GEN nf, GEN x)
     143             : {
     144     5001615 :   long tx = typ(x);
     145     5001615 :   if (tx == t_MAT || tx == t_INT) return x;
     146     4839057 :   x = nf_to_scalar_or_basis(nf, x);
     147     4838932 :   return (typ(x) == t_COL)? zk_multable(nf, x): x;
     148             : }
     149             : 
     150             : GEN
     151       21302 : nftrace(GEN nf, GEN x)
     152             : {
     153       21302 :   pari_sp av = avma;
     154       21302 :   nf = checknf(nf);
     155       21301 :   x = nf_to_scalar_or_basis(nf, x);
     156       21282 :   x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
     157       21303 :                        : gmulgu(x, nf_get_degree(nf));
     158       21304 :   return gerepileupto(av, x);
     159             : }
     160             : GEN
     161        1043 : rnfelttrace(GEN rnf, GEN x)
     162             : {
     163        1043 :   pari_sp av = avma;
     164        1043 :   checkrnf(rnf);
     165             :   /* avoid rnfabstorel special t_POL case misinterpretation */
     166        1036 :   if (typ(x) == t_POL && varn(x) == rnf_get_varn(rnf))
     167          63 :     x = gmodulo(x, rnf_get_pol(rnf));
     168        1036 :   x = rnfeltabstorel(rnf, x);
     169         721 :   x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
     170         826 :                           : gmulgu(x, rnf_get_degree(rnf));
     171         826 :   return gerepileupto(av, x);
     172             : }
     173             : 
     174             : static GEN
     175          35 : famatQ_to_famatZ(GEN fa)
     176             : {
     177          35 :   GEN E, F, Q, P = gel(fa,1);
     178          35 :   long i, j, l = lg(P);
     179          35 :   if (l == 1 || RgV_is_ZV(P)) return fa;
     180           7 :   Q = cgetg(2*l, t_COL);
     181           7 :   F = cgetg(2*l, t_COL); E = gel(fa, 2);
     182          35 :   for (i = j = 1; i < l; i++)
     183             :   {
     184          28 :     GEN p = gel(P,i);
     185          28 :     if (typ(p) == t_INT)
     186          14 :     { gel(Q, j) = p; gel(F, j) = gel(E, i); j++; }
     187             :     else
     188             :     {
     189          14 :       gel(Q, j) = gel(p,1); gel(F, j) = gel(E, i); j++;
     190          14 :       gel(Q, j) = gel(p,2); gel(F, j) = negi(gel(E, i)); j++;
     191             :     }
     192             :   }
     193           7 :   setlg(Q, j); setlg(F, j); return mkmat2(Q, F);
     194             : }
     195             : static GEN
     196          35 : famat_cba(GEN fa)
     197             : {
     198          35 :   GEN Q, F, P = gel(fa, 1), E = gel(fa, 2);
     199          35 :   long i, j, lQ, l = lg(P);
     200          35 :   if (l == 1) return fa;
     201          28 :   Q = ZV_cba(P); lQ = lg(Q); settyp(Q, t_COL);
     202          28 :   F = cgetg(lQ, t_COL);
     203          77 :   for (j = 1; j < lQ; j++)
     204             :   {
     205          49 :     GEN v = gen_0, q = gel(Q,j);
     206          49 :     if (!equali1(q))
     207         203 :       for (i = 1; i < l; i++)
     208             :       {
     209         161 :         long e = Z_pval(gel(P,i), q);
     210         161 :         if (e) v = addii(v, muliu(gel(E,i), e));
     211             :       }
     212          49 :     gel(F, j) = v;
     213             :   }
     214          28 :   return mkmat2(Q, F);
     215             : }
     216             : static long
     217          35 : famat_sign(GEN fa)
     218             : {
     219          35 :   GEN P = gel(fa,1), E = gel(fa,2);
     220          35 :   long i, l = lg(P), s = 1;
     221         126 :   for (i = 1; i < l; i++)
     222          91 :     if (signe(gel(P,i)) < 0 && mpodd(gel(E,i))) s = -s;
     223          35 :   return s;
     224             : }
     225             : static GEN
     226          35 : famat_abs(GEN fa)
     227             : {
     228          35 :   GEN Q, P = gel(fa,1);
     229             :   long i, l;
     230          35 :   Q = cgetg_copy(P, &l);
     231         126 :   for (i = 1; i < l; i++) gel(Q,i) = absi_shallow(gel(P,i));
     232          35 :   return mkmat2(Q, gel(fa,2));
     233             : }
     234             : 
     235             : /* assume nf is a genuine nf, fa a famat */
     236             : static GEN
     237          35 : famat_norm(GEN nf, GEN fa)
     238             : {
     239          35 :   pari_sp av = avma;
     240          35 :   GEN G, g = gel(fa,1);
     241             :   long i, l, s;
     242             : 
     243          35 :   G = cgetg_copy(g, &l);
     244         112 :   for (i = 1; i < l; i++) gel(G,i) = nfnorm(nf, gel(g,i));
     245          35 :   fa = mkmat2(G, gel(fa,2));
     246          35 :   fa = famatQ_to_famatZ(fa);
     247          35 :   s = famat_sign(fa);
     248          35 :   fa = famat_reduce(famat_abs(fa));
     249          35 :   fa = famat_cba(fa);
     250          35 :   g = factorback(fa);
     251          35 :   return gerepileupto(av, s < 0? gneg(g): g);
     252             : }
     253             : GEN
     254      223156 : nfnorm(GEN nf, GEN x)
     255             : {
     256      223156 :   pari_sp av = avma;
     257             :   GEN c, den;
     258             :   long n;
     259      223156 :   nf = checknf(nf);
     260      223156 :   n = nf_get_degree(nf);
     261      223156 :   if (typ(x) == t_MAT) return famat_norm(nf, x);
     262      223121 :   x = nf_to_scalar_or_basis(nf, x);
     263      223120 :   if (typ(x)!=t_COL)
     264      126910 :     return gerepileupto(av, gpowgs(x, n));
     265       96210 :   x = nf_to_scalar_or_alg(nf, Q_primitive_part(x, &c));
     266       96211 :   x = Q_remove_denom(x, &den);
     267       96211 :   x = ZX_resultant_all(nf_get_pol(nf), x, den, 0);
     268       96211 :   return gerepileupto(av, c ? gmul(x, gpowgs(c, n)): x);
     269             : }
     270             : 
     271             : static GEN
     272         119 : to_RgX(GEN P, long vx)
     273             : {
     274         119 :   return varn(P) == vx ? P: scalarpol_shallow(P, vx);
     275             : }
     276             : 
     277             : GEN
     278         462 : rnfeltnorm(GEN rnf, GEN x)
     279             : {
     280         462 :   pari_sp av = avma;
     281             :   GEN nf, pol;
     282             :   long v;
     283         462 :   checkrnf(rnf);
     284         455 :   v = rnf_get_varn(rnf);
     285             :   /* avoid rnfabstorel special t_POL case misinterpretation */
     286         455 :   if (typ(x) == t_POL && varn(x) == v) x = gmodulo(x, rnf_get_pol(rnf));
     287         455 :   x = liftpol_shallow(rnfeltabstorel(rnf, x));
     288         245 :   nf = rnf_get_nf(rnf); pol = rnf_get_pol(rnf);
     289         490 :   x = (typ(x) == t_POL)
     290         119 :     ? rnfeltdown(rnf, nfX_resultant(nf,pol,to_RgX(x,v)))
     291         245 :     : gpowgs(x, rnf_get_degree(rnf));
     292         245 :   return gerepileupto(av, x);
     293             : }
     294             : 
     295             : /* x + y in nf */
     296             : GEN
     297    23476606 : nfadd(GEN nf, GEN x, GEN y)
     298             : {
     299    23476606 :   pari_sp av = avma;
     300             :   GEN z;
     301             : 
     302    23476606 :   nf = checknf(nf);
     303    23476606 :   x = nf_to_scalar_or_basis(nf, x);
     304    23476606 :   y = nf_to_scalar_or_basis(nf, y);
     305    23476606 :   if (typ(x) != t_COL)
     306    17706842 :   { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
     307             :   else
     308     5769764 :   { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
     309    23476605 :   return gerepileupto(av, z);
     310             : }
     311             : /* x - y in nf */
     312             : GEN
     313     1815196 : nfsub(GEN nf, GEN x, GEN y)
     314             : {
     315     1815196 :   pari_sp av = avma;
     316             :   GEN z;
     317             : 
     318     1815196 :   nf = checknf(nf);
     319     1815196 :   x = nf_to_scalar_or_basis(nf, x);
     320     1815196 :   y = nf_to_scalar_or_basis(nf, y);
     321     1815196 :   if (typ(x) != t_COL)
     322     1282372 :   { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
     323             :   else
     324      532824 :   { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
     325     1815196 :   return gerepileupto(av, z);
     326             : }
     327             : 
     328             : /* product of ZC x,y in (true) nf; ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
     329             : static GEN
     330     9062107 : nfmuli_ZC(GEN nf, GEN x, GEN y)
     331             : {
     332             :   long i, j, k, N;
     333     9062107 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     334             : 
     335    44072411 :   for (k = 1; k <= N; k++)
     336             :   {
     337    35010425 :     pari_sp av = avma;
     338    35010425 :     GEN s, TABi = TAB;
     339    35010425 :     if (k == 1)
     340     9062095 :       s = mulii(gel(x,1),gel(y,1));
     341             :     else
     342    25948166 :       s = addii(mulii(gel(x,1),gel(y,k)),
     343    25948330 :                 mulii(gel(x,k),gel(y,1)));
     344   227076605 :     for (i=2; i<=N; i++)
     345             :     {
     346   192070650 :       GEN t, xi = gel(x,i);
     347   192070650 :       TABi += N;
     348   192070650 :       if (!signe(xi)) continue;
     349             : 
     350    96756994 :       t = NULL;
     351  1084254319 :       for (j=2; j<=N; j++)
     352             :       {
     353   987499214 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     354   987499214 :         if (!signe(c)) continue;
     355   291619991 :         p1 = _mulii(c, gel(y,j));
     356   291625924 :         t = t? addii(t, p1): p1;
     357             :       }
     358    96755105 :       if (t) s = addii(s, mulii(xi, t));
     359             :     }
     360    35005955 :     gel(v,k) = gerepileuptoint(av,s);
     361             :   }
     362     9061986 :   return v;
     363             : }
     364             : static int
     365    74739149 : is_famat(GEN x) { return typ(x) == t_MAT && lg(x) == 3; }
     366             : /* product of x and y in nf */
     367             : GEN
     368    36372946 : nfmul(GEN nf, GEN x, GEN y)
     369             : {
     370             :   GEN z;
     371    36372946 :   pari_sp av = avma;
     372             : 
     373    36372946 :   if (x == y) return nfsqr(nf,x);
     374             : 
     375    32276893 :   nf = checknf(nf);
     376    32276893 :   if (is_famat(x) || is_famat(y)) return famat_mul(x, y);
     377    32276584 :   x = nf_to_scalar_or_basis(nf, x);
     378    32276585 :   y = nf_to_scalar_or_basis(nf, y);
     379    32276580 :   if (typ(x) != t_COL)
     380             :   {
     381    21843941 :     if (isintzero(x)) return gen_0;
     382    15772813 :     z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
     383             :   else
     384             :   {
     385    10432639 :     if (typ(y) != t_COL)
     386             :     {
     387     4547657 :       if (isintzero(y)) return gen_0;
     388     1613539 :       z = RgC_Rg_mul(x, y);
     389             :     }
     390             :     else
     391             :     {
     392             :       GEN dx, dy;
     393     5884982 :       x = Q_remove_denom(x, &dx);
     394     5884982 :       y = Q_remove_denom(y, &dy);
     395     5884983 :       z = nfmuli_ZC(nf,x,y);
     396     5884981 :       dx = mul_denom(dx,dy);
     397     5884981 :       if (dx) z = ZC_Z_div(z, dx);
     398             :     }
     399             :   }
     400    23271331 :   return gerepileupto(av, z);
     401             : }
     402             : /* square of ZC x in nf */
     403             : static GEN
     404     7129215 : nfsqri_ZC(GEN nf, GEN x)
     405             : {
     406             :   long i, j, k, N;
     407     7129215 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     408             : 
     409    39014409 :   for (k = 1; k <= N; k++)
     410             :   {
     411    31885248 :     pari_sp av = avma;
     412    31885248 :     GEN s, TABi = TAB;
     413    31885248 :     if (k == 1)
     414     7129394 :       s = sqri(gel(x,1));
     415             :     else
     416    24755854 :       s = shifti(mulii(gel(x,1),gel(x,k)), 1);
     417   253838524 :     for (i=2; i<=N; i++)
     418             :     {
     419   221973486 :       GEN p1, c, t, xi = gel(x,i);
     420   221973486 :       TABi += N;
     421   221973486 :       if (!signe(xi)) continue;
     422             : 
     423    79970232 :       c = gcoeff(TABi, k, i);
     424    79970232 :       t = signe(c)? _mulii(c,xi): NULL;
     425   676102514 :       for (j=i+1; j<=N; j++)
     426             :       {
     427   596131344 :         c = gcoeff(TABi, k, j);
     428   596131344 :         if (!signe(c)) continue;
     429   231963096 :         p1 = _mulii(c, shifti(gel(x,j),1));
     430   231970592 :         t = t? addii(t, p1): p1;
     431             :       }
     432    79971170 :       if (t) s = addii(s, mulii(xi, t));
     433             :     }
     434    31865038 :     gel(v,k) = gerepileuptoint(av,s);
     435             :   }
     436     7129161 :   return v;
     437             : }
     438             : /* square of x in nf */
     439             : GEN
     440     8914784 : nfsqr(GEN nf, GEN x)
     441             : {
     442     8914784 :   pari_sp av = avma;
     443             :   GEN z;
     444             : 
     445     8914784 :   nf = checknf(nf);
     446     8914789 :   if (is_famat(x)) return famat_sqr(x);
     447     8914792 :   x = nf_to_scalar_or_basis(nf, x);
     448     8914793 :   if (typ(x) != t_COL) z = gsqr(x);
     449             :   else
     450             :   {
     451             :     GEN dx;
     452     2632090 :     x = Q_remove_denom(x, &dx);
     453     2632091 :     z = nfsqri_ZC(nf,x);
     454     2632089 :     if (dx) z = RgC_Rg_div(z, sqri(dx));
     455             :   }
     456     8914792 :   return gerepileupto(av, z);
     457             : }
     458             : 
     459             : /* x a ZC, v a t_COL of ZC/Z */
     460             : GEN
     461      205721 : zkC_multable_mul(GEN v, GEN x)
     462             : {
     463      205721 :   long i, l = lg(v);
     464      205721 :   GEN y = cgetg(l, t_COL);
     465      800269 :   for (i = 1; i < l; i++)
     466             :   {
     467      594548 :     GEN c = gel(v,i);
     468      594548 :     if (typ(c)!=t_COL) {
     469           0 :       if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
     470             :     } else {
     471      594548 :       c = ZM_ZC_mul(x,c);
     472      594548 :       if (ZV_isscalar(c)) c = gel(c,1);
     473             :     }
     474      594548 :     gel(y,i) = c;
     475             :   }
     476      205721 :   return y;
     477             : }
     478             : 
     479             : GEN
     480       57227 : nfC_multable_mul(GEN v, GEN x)
     481             : {
     482       57227 :   long i, l = lg(v);
     483       57227 :   GEN y = cgetg(l, t_COL);
     484      385363 :   for (i = 1; i < l; i++)
     485             :   {
     486      328136 :     GEN c = gel(v,i);
     487      328136 :     if (typ(c)!=t_COL) {
     488      273526 :       if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
     489             :     } else {
     490       54610 :       c = RgM_RgC_mul(x,c);
     491       54610 :       if (QV_isscalar(c)) c = gel(c,1);
     492             :     }
     493      328136 :     gel(y,i) = c;
     494             :   }
     495       57227 :   return y;
     496             : }
     497             : 
     498             : GEN
     499      200022 : nfC_nf_mul(GEN nf, GEN v, GEN x)
     500             : {
     501             :   long tx;
     502             :   GEN y;
     503             : 
     504      200022 :   x = nf_to_scalar_or_basis(nf, x);
     505      200022 :   tx = typ(x);
     506      200022 :   if (tx != t_COL)
     507             :   {
     508             :     long l, i;
     509      151425 :     if (tx == t_INT)
     510             :     {
     511      142150 :       long s = signe(x);
     512      142150 :       if (!s) return zerocol(lg(v)-1);
     513      134676 :       if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
     514             :     }
     515       49098 :     l = lg(v); y = cgetg(l, t_COL);
     516      350483 :     for (i=1; i < l; i++)
     517             :     {
     518      301385 :       GEN c = gel(v,i);
     519      301385 :       if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
     520      301385 :       gel(y,i) = c;
     521             :     }
     522       49098 :     return y;
     523             :   }
     524             :   else
     525             :   {
     526             :     GEN dx;
     527       48597 :     x = zk_multable(nf, Q_remove_denom(x,&dx));
     528       48597 :     y = nfC_multable_mul(v, x);
     529       48597 :     return dx? RgC_Rg_div(y, dx): y;
     530             :   }
     531             : }
     532             : static GEN
     533       11213 : mulbytab(GEN M, GEN c)
     534       11213 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
     535             : GEN
     536        2751 : tablemulvec(GEN M, GEN x, GEN v)
     537             : {
     538             :   long l, i;
     539             :   GEN y;
     540             : 
     541        2751 :   if (typ(x) == t_COL && RgV_isscalar(x))
     542             :   {
     543           0 :     x = gel(x,1);
     544           0 :     return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
     545             :   }
     546        2751 :   x = multable(M, x); /* multiplication table by x */
     547        2751 :   y = cgetg_copy(v, &l);
     548        2751 :   if (typ(v) == t_POL)
     549             :   {
     550        2751 :     y[1] = v[1];
     551       13964 :     for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     552        2751 :     y = normalizepol(y);
     553             :   }
     554             :   else
     555             :   {
     556           0 :     for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     557             :   }
     558        2751 :   return y;
     559             : }
     560             : 
     561             : GEN
     562     1261705 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
     563             : GEN
     564     1580945 : zkmultable_inv(GEN mx) { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
     565             : /* nf a true nf, x a ZC */
     566             : GEN
     567      319244 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
     568             : 
     569             : /* inverse of x in nf */
     570             : GEN
     571      240289 : nfinv(GEN nf, GEN x)
     572             : {
     573      240289 :   pari_sp av = avma;
     574             :   GEN z;
     575             : 
     576      240289 :   nf = checknf(nf);
     577      240289 :   if (is_famat(x)) return famat_inv(x);
     578      240289 :   x = nf_to_scalar_or_basis(nf, x);
     579      240289 :   if (typ(x) == t_COL)
     580             :   {
     581             :     GEN d;
     582      190865 :     x = Q_remove_denom(x, &d);
     583      190865 :     z = zk_inv(nf, x);
     584      190865 :     if (d) z = RgC_Rg_mul(z, d);
     585             :   }
     586             :   else
     587       49424 :     z = ginv(x);
     588      240289 :   return gerepileupto(av, z);
     589             : }
     590             : 
     591             : /* quotient of x and y in nf */
     592             : GEN
     593       36319 : nfdiv(GEN nf, GEN x, GEN y)
     594             : {
     595       36319 :   pari_sp av = avma;
     596             :   GEN z;
     597             : 
     598       36319 :   nf = checknf(nf);
     599       36319 :   if (is_famat(x) || is_famat(y)) return famat_div(x,y);
     600       36228 :   y = nf_to_scalar_or_basis(nf, y);
     601       36228 :   if (typ(y) != t_COL)
     602             :   {
     603       22099 :     x = nf_to_scalar_or_basis(nf, x);
     604       22099 :     z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
     605             :   }
     606             :   else
     607             :   {
     608             :     GEN d;
     609       14129 :     y = Q_remove_denom(y, &d);
     610       14128 :     z = nfmul(nf, x, zk_inv(nf,y));
     611       14129 :     if (d) z = typ(z) == t_COL? RgC_Rg_mul(z, d): gmul(z, d);
     612             :   }
     613       36228 :   return gerepileupto(av, z);
     614             : }
     615             : 
     616             : /* product of INTEGERS (t_INT or ZC) x and y in (true) nf */
     617             : GEN
     618     4549369 : nfmuli(GEN nf, GEN x, GEN y)
     619             : {
     620     4549369 :   if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
     621     3410123 :   if (typ(y) == t_INT) return ZC_Z_mul(x, y);
     622     3177099 :   return nfmuli_ZC(nf, x, y);
     623             : }
     624             : GEN
     625     4497170 : nfsqri(GEN nf, GEN x)
     626     4497170 : { return (typ(x) == t_INT)? sqri(x): nfsqri_ZC(nf, x); }
     627             : 
     628             : /* both x and y are RgV */
     629             : GEN
     630           0 : tablemul(GEN TAB, GEN x, GEN y)
     631             : {
     632             :   long i, j, k, N;
     633             :   GEN s, v;
     634           0 :   if (typ(x) != t_COL) return gmul(x, y);
     635           0 :   if (typ(y) != t_COL) return gmul(y, x);
     636           0 :   N = lg(x)-1;
     637           0 :   v = cgetg(N+1,t_COL);
     638           0 :   for (k=1; k<=N; k++)
     639             :   {
     640           0 :     pari_sp av = avma;
     641           0 :     GEN TABi = TAB;
     642           0 :     if (k == 1)
     643           0 :       s = gmul(gel(x,1),gel(y,1));
     644             :     else
     645           0 :       s = gadd(gmul(gel(x,1),gel(y,k)),
     646           0 :                gmul(gel(x,k),gel(y,1)));
     647           0 :     for (i=2; i<=N; i++)
     648             :     {
     649           0 :       GEN t, xi = gel(x,i);
     650           0 :       TABi += N;
     651           0 :       if (gequal0(xi)) continue;
     652             : 
     653           0 :       t = NULL;
     654           0 :       for (j=2; j<=N; j++)
     655             :       {
     656           0 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     657           0 :         if (gequal0(c)) continue;
     658           0 :         p1 = gmul(c, gel(y,j));
     659           0 :         t = t? gadd(t, p1): p1;
     660             :       }
     661           0 :       if (t) s = gadd(s, gmul(xi, t));
     662             :     }
     663           0 :     gel(v,k) = gerepileupto(av,s);
     664             :   }
     665           0 :   return v;
     666             : }
     667             : GEN
     668       49090 : tablesqr(GEN TAB, GEN x)
     669             : {
     670             :   long i, j, k, N;
     671             :   GEN s, v;
     672             : 
     673       49090 :   if (typ(x) != t_COL) return gsqr(x);
     674       49090 :   N = lg(x)-1;
     675       49090 :   v = cgetg(N+1,t_COL);
     676             : 
     677      349430 :   for (k=1; k<=N; k++)
     678             :   {
     679      300340 :     pari_sp av = avma;
     680      300340 :     GEN TABi = TAB;
     681      300340 :     if (k == 1)
     682       49090 :       s = gsqr(gel(x,1));
     683             :     else
     684      251250 :       s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
     685     1909862 :     for (i=2; i<=N; i++)
     686             :     {
     687     1609522 :       GEN p1, c, t, xi = gel(x,i);
     688     1609522 :       TABi += N;
     689     1609522 :       if (gequal0(xi)) continue;
     690             : 
     691      419846 :       c = gcoeff(TABi, k, i);
     692      419846 :       t = !gequal0(c)? gmul(c,xi): NULL;
     693     1676969 :       for (j=i+1; j<=N; j++)
     694             :       {
     695     1257123 :         c = gcoeff(TABi, k, j);
     696     1257123 :         if (gequal0(c)) continue;
     697      646443 :         p1 = gmul(gmul2n(c,1), gel(x,j));
     698      646443 :         t = t? gadd(t, p1): p1;
     699             :       }
     700      419846 :       if (t) s = gadd(s, gmul(xi, t));
     701             :     }
     702      300340 :     gel(v,k) = gerepileupto(av,s);
     703             :   }
     704       49090 :   return v;
     705             : }
     706             : 
     707             : static GEN
     708      356427 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
     709             : static GEN
     710      987162 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
     711             : 
     712             : /* Compute z^n in nf, left-shift binary powering */
     713             : GEN
     714      942997 : nfpow(GEN nf, GEN z, GEN n)
     715             : {
     716      942997 :   pari_sp av = avma;
     717             :   long s;
     718             :   GEN x, cx;
     719             : 
     720      942997 :   if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
     721      942997 :   nf = checknf(nf);
     722      942995 :   s = signe(n); if (!s) return gen_1;
     723      942995 :   if (is_famat(z)) return famat_pow(z, n);
     724      882354 :   x = nf_to_scalar_or_basis(nf, z);
     725      882357 :   if (typ(x) != t_COL) return powgi(x,n);
     726      762524 :   if (s < 0)
     727             :   { /* simplified nfinv */
     728             :     GEN d;
     729       45754 :     x = Q_remove_denom(x, &d);
     730       45754 :     x = zk_inv(nf, x);
     731       45754 :     x = primitive_part(x, &cx);
     732       45753 :     cx = mul_content(cx, d);
     733       45753 :     n = negi(n);
     734             :   }
     735             :   else
     736      716770 :     x = primitive_part(x, &cx);
     737      762499 :   x = gen_pow_i(x, n, (void*)nf, _sqr, _mul);
     738      762514 :   if (cx)
     739       46791 :     x = gerepileupto(av, gmul(x, powgi(cx, n)));
     740             :   else
     741      715723 :     x = gerepilecopy(av, x);
     742      762533 :   return x;
     743             : }
     744             : /* Compute z^n in nf, left-shift binary powering */
     745             : GEN
     746      354716 : nfpow_u(GEN nf, GEN z, ulong n)
     747             : {
     748      354716 :   pari_sp av = avma;
     749             :   GEN x, cx;
     750             : 
     751      354716 :   if (!n) return gen_1;
     752      354716 :   x = nf_to_scalar_or_basis(nf, z);
     753      354716 :   if (typ(x) != t_COL) return gpowgs(x,n);
     754      318453 :   x = primitive_part(x, &cx);
     755      318453 :   x = gen_powu_i(x, n, (void*)nf, _sqr, _mul);
     756      318451 :   if (cx)
     757             :   {
     758      114520 :     x = gmul(x, powgi(cx, utoipos(n)));
     759      114520 :     return gerepileupto(av,x);
     760             :   }
     761      203931 :   return gerepilecopy(av, x);
     762             : }
     763             : 
     764             : long
     765        1099 : nfissquare(GEN nf, GEN z, GEN *px)
     766             : {
     767        1099 :   pari_sp av = avma;
     768        1099 :   long v = fetch_var_higher();
     769             :   GEN R;
     770        1099 :   nf = checknf(nf);
     771        1099 :   if (nf_get_degree(nf) == 1)
     772             :   {
     773         189 :     z = algtobasis(nf, z);
     774         189 :     if (!issquareall(gel(z,1), px)) return gc_long(av, 0);
     775          21 :     if (px) *px = gerepileupto(av, *px); else set_avma(av);
     776          21 :     return 1;
     777             :   }
     778         910 :   z = nf_to_scalar_or_alg(nf, z);
     779         910 :   R = nfroots(nf, deg2pol_shallow(gen_m1, gen_0, z, v));
     780         910 :   delete_var(); if (lg(R) == 1) return gc_long(av, 0);
     781         560 :   if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
     782          14 :   else set_avma(av);
     783         560 :   return 1;
     784             : }
     785             : 
     786             : long
     787        7720 : nfispower(GEN nf, GEN z, long n, GEN *px)
     788             : {
     789        7720 :   pari_sp av = avma;
     790        7720 :   long v = fetch_var_higher();
     791             :   GEN R;
     792        7720 :   nf = checknf(nf);
     793        7720 :   if (nf_get_degree(nf) == 1)
     794             :   {
     795         329 :     z = algtobasis(nf, z);
     796         329 :     if (!ispower(gel(z,1), stoi(n), px)) return gc_long(av, 0);
     797         147 :     if (px) *px = gerepileupto(av, *px); else set_avma(av);
     798         147 :     return 1;
     799             :   }
     800        7391 :   if (n <= 0)
     801           0 :     pari_err_DOMAIN("nfeltispower","exponent","<=",gen_0,stoi(n));
     802        7391 :   z = nf_to_scalar_or_alg(nf, z);
     803        7391 :   if (n==1)
     804             :   {
     805           0 :     if (px) *px = gerepilecopy(av, z);
     806           0 :     return 1;
     807             :   }
     808        7391 :   R = nfroots(nf, gsub(pol_xn(n, v), z));
     809        7391 :   delete_var(); if (lg(R) == 1) return gc_long(av, 0);
     810        3157 :   if (px) *px = gerepilecopy(av, nf_to_scalar_or_basis(nf, gel(R,1)));
     811        3143 :   else set_avma(av);
     812        3157 :   return 1;
     813             : }
     814             : 
     815             : static GEN
     816          56 : idmulred(void *nf, GEN x, GEN y) { return idealmulred((GEN) nf, x, y); }
     817             : static GEN
     818         413 : idpowred(void *nf, GEN x, GEN n) { return idealpowred((GEN) nf, x, n); }
     819             : static GEN
     820       72020 : idmul(void *nf, GEN x, GEN y) { return idealmul((GEN) nf, x, y); }
     821             : static GEN
     822       87971 : idpow(void *nf, GEN x, GEN n) { return idealpow((GEN) nf, x, n); }
     823             : GEN
     824       86368 : idealfactorback(GEN nf, GEN L, GEN e, long red)
     825             : {
     826       86368 :   nf = checknf(nf);
     827       86368 :   if (red) return gen_factorback(L, e, (void*)nf, &idmulred, &idpowred, NULL);
     828       86011 :   if (!e && typ(L) == t_MAT && lg(L) == 3) { e = gel(L,2); L = gel(L,1); }
     829       86011 :   if (is_vec_t(typ(L)) && RgV_is_prV(L))
     830             :   { /* don't use gen_factorback since *= pr^v can be done more efficiently */
     831       65377 :     pari_sp av = avma;
     832       65377 :     long i, l = lg(L);
     833             :     GEN a;
     834       65377 :     if (!e) e = const_vec(l-1, gen_1);
     835       62521 :     else switch(typ(e))
     836             :     {
     837        7768 :       case t_VECSMALL: e = zv_to_ZV(e); break;
     838       54753 :       case t_VEC: case t_COL:
     839       54753 :         if (!RgV_is_ZV(e))
     840           0 :           pari_err_TYPE("factorback [not an exponent vector]", e);
     841       54753 :         break;
     842           0 :       default: pari_err_TYPE("idealfactorback", e);
     843             :     }
     844       65377 :     if (l != lg(e))
     845           0 :       pari_err_TYPE("factorback [not an exponent vector]", e);
     846       65377 :     if (l == 1 || ZV_equal0(e)) return gc_const(av, gen_1);
     847       23714 :     a = idealpow(nf, gel(L,1), gel(e,1));
     848      252093 :     for (i = 2; i < l; i++)
     849      228379 :       if (signe(gel(e,i))) a = idealmulpowprime(nf, a, gel(L,i), gel(e,i));
     850       23714 :     return gerepileupto(av, a);
     851             :   }
     852       20634 :   return gen_factorback(L, e, (void*)nf, &idmul, &idpow, NULL);
     853             : }
     854             : static GEN
     855      327935 : eltmul(void *nf, GEN x, GEN y) { return nfmul((GEN) nf, x, y); }
     856             : static GEN
     857      465340 : eltpow(void *nf, GEN x, GEN n) { return nfpow((GEN) nf, x, n); }
     858             : GEN
     859      265468 : nffactorback(GEN nf, GEN L, GEN e)
     860      265468 : { return gen_factorback(L, e, (void*)checknf(nf), &eltmul, &eltpow, NULL); }
     861             : 
     862             : static GEN
     863     3099304 : _nf_red(void *E, GEN x) { (void)E; return gcopy(x); }
     864             : 
     865             : static GEN
     866    12672932 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
     867             : 
     868             : static GEN
     869      751655 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
     870             : 
     871             : static GEN
     872    15218661 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
     873             : 
     874             : static GEN
     875       53959 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
     876             : 
     877             : static GEN
     878       11128 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
     879             : 
     880             : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
     881             :                                         _nf_inv,&gequal0,_nf_s };
     882             : 
     883      227960 : const struct bb_field *get_nf_field(void **E, GEN nf)
     884      227960 : { *E = (void*)nf; return &nf_field; }
     885             : 
     886             : GEN
     887          14 : nfM_det(GEN nf, GEN M)
     888             : {
     889             :   void *E;
     890          14 :   const struct bb_field *S = get_nf_field(&E, nf);
     891          14 :   return gen_det(M, E, S);
     892             : }
     893             : GEN
     894       11114 : nfM_inv(GEN nf, GEN M)
     895             : {
     896             :   void *E;
     897       11114 :   const struct bb_field *S = get_nf_field(&E, nf);
     898       11114 :   return gen_Gauss(M, matid(lg(M)-1), E, S);
     899             : }
     900             : 
     901             : GEN
     902           0 : nfM_ker(GEN nf, GEN M)
     903             : {
     904             :    void *E;
     905           0 :    const struct bb_field *S = get_nf_field(&E, nf);
     906           0 :    return gen_ker(M, 0, E, S);
     907             : }
     908             : 
     909             : GEN
     910       10610 : nfM_mul(GEN nf, GEN A, GEN B)
     911             : {
     912             :   void *E;
     913       10610 :   const struct bb_field *S = get_nf_field(&E, nf);
     914       10610 :   return gen_matmul(A, B, E, S);
     915             : }
     916             : GEN
     917      206222 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
     918             : {
     919             :   void *E;
     920      206222 :   const struct bb_field *S = get_nf_field(&E, nf);
     921      206222 :   return gen_matcolmul(A, B, E, S);
     922             : }
     923             : 
     924             : /* valuation of integral x (ZV), with resp. to prime ideal pr */
     925             : long
     926    24032287 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
     927             : {
     928    24032287 :   pari_sp av = avma;
     929             :   long i, v, l;
     930    24032287 :   GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
     931             : 
     932             :   /* p inert */
     933    24032305 :   if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
     934    23026287 :   y = cgetg_copy(x, &l); /* will hold the new x */
     935    23026658 :   x = leafcopy(x);
     936    37197715 :   for(v=0;; v++)
     937             :   {
     938   143064916 :     for (i=1; i<l; i++)
     939             :     { /* is (x.b)[i] divisible by p ? */
     940   128888171 :       gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
     941   128891752 :       if (r != gen_0) { if (newx) *newx = x; return v; }
     942             :     }
     943    14176745 :     swap(x, y);
     944    14176745 :     if (!newx && (v & 0xf) == 0xf) v += pr_get_e(pr) * ZV_pvalrem(x, p, &x);
     945    14176745 :     if (gc_needed(av,1))
     946             :     {
     947           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZC_nfvalrem, v >= %ld", v);
     948           0 :       gerepileall(av, 2, &x, &y);
     949             :     }
     950             :   }
     951             : }
     952             : long
     953    19756089 : ZC_nfval(GEN x, GEN P)
     954    19756089 : { return ZC_nfvalrem(x, P, NULL); }
     955             : 
     956             : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
     957             : int
     958     1250184 : ZC_prdvd(GEN x, GEN P)
     959             : {
     960     1250184 :   pari_sp av = avma;
     961             :   long i, l;
     962     1250184 :   GEN p = pr_get_p(P), mul = pr_get_tau(P);
     963     1250195 :   if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
     964     1249649 :   l = lg(x);
     965     5063761 :   for (i=1; i<l; i++)
     966     4546464 :     if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) return gc_bool(av,0);
     967      517297 :   return gc_bool(av,1);
     968             : }
     969             : 
     970             : int
     971         357 : pr_equal(GEN P, GEN Q)
     972             : {
     973         357 :   GEN gQ, p = pr_get_p(P);
     974         357 :   long e = pr_get_e(P), f = pr_get_f(P), n;
     975         357 :   if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
     976         336 :     return 0;
     977          21 :   gQ = pr_get_gen(Q); n = lg(gQ)-1;
     978          21 :   if (2*e*f > n) return 1; /* room for only one such pr */
     979          14 :   return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
     980             : }
     981             : 
     982             : GEN
     983      420721 : famat_nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     984             : {
     985      420721 :   pari_sp av = avma;
     986      420721 :   GEN P = gel(x,1), E = gel(x,2), V = gen_0, y = NULL;
     987      420721 :   long l = lg(P), simplify = 0, i;
     988      420721 :   if (py) { *py = gen_1; y = cgetg(l, t_COL); }
     989             : 
     990     2259113 :   for (i = 1; i < l; i++)
     991             :   {
     992     1838392 :     GEN e = gel(E,i);
     993             :     long v;
     994     1838392 :     if (!signe(e))
     995             :     {
     996           7 :       if (py) gel(y,i) = gen_1;
     997           7 :       simplify = 1; continue;
     998             :     }
     999     1838385 :     v = nfvalrem(nf, gel(P,i), pr, py? &gel(y,i): NULL);
    1000     1838385 :     if (v == LONG_MAX) { set_avma(av); if (py) *py = gen_0; return mkoo(); }
    1001     1838385 :     V = addmulii(V, stoi(v), e);
    1002             :   }
    1003      420721 :   if (!py) V = gerepileuptoint(av, V);
    1004             :   else
    1005             :   {
    1006          42 :     y = mkmat2(y, gel(x,2));
    1007          42 :     if (simplify) y = famat_remove_trivial(y);
    1008          42 :     gerepileall(av, 2, &V, &y); *py = y;
    1009             :   }
    1010      420721 :   return V;
    1011             : }
    1012             : long
    1013     5633058 : nfval(GEN nf, GEN x, GEN pr)
    1014             : {
    1015     5633058 :   pari_sp av = avma;
    1016             :   long w, e;
    1017             :   GEN cx, p;
    1018             : 
    1019     5633058 :   if (gequal0(x)) return LONG_MAX;
    1020     5619601 :   nf = checknf(nf);
    1021     5619601 :   checkprid(pr);
    1022     5619589 :   p = pr_get_p(pr);
    1023     5619587 :   e = pr_get_e(pr);
    1024     5619583 :   x = nf_to_scalar_or_basis(nf, x);
    1025     5619488 :   if (typ(x) != t_COL) return e*Q_pval(x,p);
    1026     2380959 :   x = Q_primitive_part(x, &cx);
    1027     2381023 :   w = ZC_nfval(x,pr);
    1028     2380946 :   if (cx) w += e*Q_pval(cx,p);
    1029     2380954 :   return gc_long(av,w);
    1030             : }
    1031             : 
    1032             : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
    1033             : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
    1034             : static GEN
    1035      973406 : powp(GEN nf, GEN pr, long v)
    1036             : {
    1037             :   GEN b, z;
    1038             :   long e;
    1039      973406 :   if (!v) return gen_1;
    1040      446810 :   b = pr_get_tau(pr);
    1041      446810 :   if (typ(b) == t_INT) return gen_1;
    1042      131320 :   e = pr_get_e(pr);
    1043      131320 :   z = gel(b,1);
    1044      131320 :   if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
    1045      131320 :   if (v < 0) { v = -v; z = nfinv(nf, z); }
    1046      131320 :   if (v != 1) z = nfpow_u(nf, z, v);
    1047      131320 :   return z;
    1048             : }
    1049             : long
    1050     3662645 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
    1051             : {
    1052     3662645 :   pari_sp av = avma;
    1053             :   long w, e;
    1054             :   GEN cx, p, t;
    1055             : 
    1056     3662645 :   if (!py) return nfval(nf,x,pr);
    1057     1810869 :   if (gequal0(x)) { *py = gen_0; return LONG_MAX; }
    1058     1810812 :   nf = checknf(nf);
    1059     1810812 :   checkprid(pr);
    1060     1810812 :   p = pr_get_p(pr);
    1061     1810812 :   e = pr_get_e(pr);
    1062     1810812 :   x = nf_to_scalar_or_basis(nf, x);
    1063     1810812 :   if (typ(x) != t_COL) {
    1064      557851 :     w = Q_pvalrem(x,p, py);
    1065      557851 :     if (!w) { *py = gerepilecopy(av, x); return 0; }
    1066      349272 :     *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
    1067      349272 :     return e*w;
    1068             :   }
    1069     1252961 :   x = Q_primitive_part(x, &cx);
    1070     1252959 :   w = ZC_nfvalrem(x,pr, py);
    1071     1252948 :   if (cx)
    1072             :   {
    1073      624134 :     long v = Q_pvalrem(cx,p, &t);
    1074      624134 :     *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
    1075      624134 :     *py = gerepileupto(av, *py);
    1076      624134 :     w += e*v;
    1077             :   }
    1078             :   else
    1079      628814 :     *py = gerepilecopy(av, *py);
    1080     1252965 :   return w;
    1081             : }
    1082             : GEN
    1083       15015 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
    1084             : {
    1085             :   long v;
    1086       15015 :   if (is_famat(x)) return famat_nfvalrem(nf, x, pr, py);
    1087       15008 :   v = nfvalrem(nf,x,pr,py);
    1088       15008 :   return v == LONG_MAX? mkoo(): stoi(v);
    1089             : }
    1090             : 
    1091             : /* true nf */
    1092             : GEN
    1093      335621 : coltoalg(GEN nf, GEN x)
    1094             : {
    1095      335621 :   return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
    1096             : }
    1097             : 
    1098             : GEN
    1099      406031 : basistoalg(GEN nf, GEN x)
    1100             : {
    1101             :   GEN T;
    1102             : 
    1103      406031 :   nf = checknf(nf);
    1104      406031 :   switch(typ(x))
    1105             :   {
    1106      329356 :     case t_COL: {
    1107      329356 :       pari_sp av = avma;
    1108      329356 :       return gerepilecopy(av, coltoalg(nf, x));
    1109             :     }
    1110       40768 :     case t_POLMOD:
    1111       40768 :       T = nf_get_pol(nf);
    1112       40768 :       if (!RgX_equal_var(T,gel(x,1)))
    1113           0 :         pari_err_MODULUS("basistoalg", T,gel(x,1));
    1114       40768 :       return gcopy(x);
    1115        6307 :     case t_POL:
    1116        6307 :       T = nf_get_pol(nf);
    1117        6307 :       if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
    1118        6300 :       retmkpolmod(RgX_rem(x, T), ZX_copy(T));
    1119       29600 :     case t_INT:
    1120             :     case t_FRAC:
    1121       29600 :       T = nf_get_pol(nf);
    1122       29600 :       retmkpolmod(gcopy(x), ZX_copy(T));
    1123           0 :     default:
    1124           0 :       pari_err_TYPE("basistoalg",x);
    1125             :       return NULL; /* LCOV_EXCL_LINE */
    1126             :   }
    1127             : }
    1128             : 
    1129             : /* true nf, x a t_POL */
    1130             : static GEN
    1131     4590310 : pol_to_scalar_or_basis(GEN nf, GEN x)
    1132             : {
    1133     4590310 :   GEN T = nf_get_pol(nf);
    1134     4590307 :   long l = lg(x);
    1135     4590307 :   if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
    1136     4590202 :   if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
    1137     4590202 :   if (l == 2) return gen_0;
    1138     3578247 :   if (l == 3)
    1139             :   {
    1140      839252 :     x = gel(x,2);
    1141      839252 :     if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
    1142      839245 :     return x;
    1143             :   }
    1144     2738995 :   return poltobasis(nf,x);
    1145             : }
    1146             : /* Assume nf is a genuine nf. */
    1147             : GEN
    1148   162241448 : nf_to_scalar_or_basis(GEN nf, GEN x)
    1149             : {
    1150   162241448 :   switch(typ(x))
    1151             :   {
    1152    97667902 :     case t_INT: case t_FRAC:
    1153    97667902 :       return x;
    1154      565059 :     case t_POLMOD:
    1155      565059 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
    1156      564922 :       switch(typ(x))
    1157             :       {
    1158       85848 :         case t_INT: case t_FRAC: return x;
    1159      479074 :         case t_POL: return pol_to_scalar_or_basis(nf,x);
    1160             :       }
    1161           0 :       break;
    1162     4111236 :     case t_POL: return pol_to_scalar_or_basis(nf,x);
    1163    59901173 :     case t_COL:
    1164    59901173 :       if (lg(x)-1 != nf_get_degree(nf)) break;
    1165    59900911 :       return QV_isscalar(x)? gel(x,1): x;
    1166             :   }
    1167          96 :   pari_err_TYPE("nf_to_scalar_or_basis",x);
    1168             :   return NULL; /* LCOV_EXCL_LINE */
    1169             : }
    1170             : /* Let x be a polynomial with coefficients in Q or nf. Return the same
    1171             :  * polynomial with coefficients expressed as vectors (on the integral basis).
    1172             :  * No consistency checks, not memory-clean. */
    1173             : GEN
    1174       29231 : RgX_to_nfX(GEN nf, GEN x)
    1175             : {
    1176             :   long i, l;
    1177       29231 :   GEN y = cgetg_copy(x, &l); y[1] = x[1];
    1178      237541 :   for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
    1179       29231 :   return y;
    1180             : }
    1181             : 
    1182             : /* Assume nf is a genuine nf. */
    1183             : GEN
    1184     4825357 : nf_to_scalar_or_alg(GEN nf, GEN x)
    1185             : {
    1186     4825357 :   switch(typ(x))
    1187             :   {
    1188       85258 :     case t_INT: case t_FRAC:
    1189       85258 :       return x;
    1190         427 :     case t_POLMOD:
    1191         427 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
    1192         427 :       if (typ(x) != t_POL) return x;
    1193             :       /* fall through */
    1194             :     case t_POL:
    1195             :     {
    1196        5334 :       GEN T = nf_get_pol(nf);
    1197        5334 :       long l = lg(x);
    1198        5334 :       if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
    1199        5334 :       if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
    1200        5334 :       if (l == 2) return gen_0;
    1201        5334 :       if (l == 3) return gel(x,2);
    1202        3794 :       return x;
    1203             :     }
    1204     4734772 :     case t_COL:
    1205             :     {
    1206             :       GEN dx;
    1207     4734772 :       if (lg(x)-1 != nf_get_degree(nf)) break;
    1208     9374914 :       if (QV_isscalar(x)) return gel(x,1);
    1209     4639945 :       x = Q_remove_denom(x, &dx);
    1210     4639986 :       x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
    1211     4640177 :       dx = mul_denom(dx, nf_get_zkden(nf));
    1212     4640147 :       return gdiv(x,dx);
    1213             :     }
    1214             :   }
    1215          54 :   pari_err_TYPE("nf_to_scalar_or_alg",x);
    1216             :   return NULL; /* LCOV_EXCL_LINE */
    1217             : }
    1218             : 
    1219             : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
    1220             : GEN
    1221        1365 : RgM_RgX_mul(GEN A, GEN x)
    1222             : {
    1223        1365 :   long i, l = lg(x)-1;
    1224             :   GEN z;
    1225        1365 :   if (l == 1) return zerocol(nbrows(A));
    1226        1351 :   z = gmul(gel(x,2), gel(A,1));
    1227        2555 :   for (i = 2; i < l; i++)
    1228        1204 :     if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
    1229        1351 :   return z;
    1230             : }
    1231             : GEN
    1232    10366446 : ZM_ZX_mul(GEN A, GEN x)
    1233             : {
    1234    10366446 :   long i, l = lg(x)-1;
    1235             :   GEN z;
    1236    10366446 :   if (l == 1) return zerocol(nbrows(A));
    1237    10365312 :   z = ZC_Z_mul(gel(A,1), gel(x,2));
    1238    32341633 :   for (i = 2; i < l ; i++)
    1239    21978906 :     if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
    1240    10362727 :   return z;
    1241             : }
    1242             : /* x a t_POL, nf a genuine nf. No garbage collecting. No check.  */
    1243             : GEN
    1244     9766802 : poltobasis(GEN nf, GEN x)
    1245             : {
    1246     9766802 :   GEN d, T = nf_get_pol(nf);
    1247     9766732 :   if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
    1248     9766599 :   if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1249     9766617 :   x = Q_remove_denom(x, &d);
    1250     9766767 :   if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
    1251     9766686 :   x = ZM_ZX_mul(nf_get_invzk(nf), x);
    1252     9764733 :   if (d) x = RgC_Rg_div(x, d);
    1253     9764799 :   return x;
    1254             : }
    1255             : 
    1256             : GEN
    1257      952857 : algtobasis(GEN nf, GEN x)
    1258             : {
    1259             :   pari_sp av;
    1260             : 
    1261      952857 :   nf = checknf(nf);
    1262      952857 :   switch(typ(x))
    1263             :   {
    1264      140547 :     case t_POLMOD:
    1265      140547 :       if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
    1266           7 :         pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
    1267      140540 :       x = gel(x,2);
    1268      140540 :       switch(typ(x))
    1269             :       {
    1270       11340 :         case t_INT:
    1271       11340 :         case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1272      129200 :         case t_POL:
    1273      129200 :           av = avma;
    1274      129200 :           return gerepileupto(av,poltobasis(nf,x));
    1275             :       }
    1276           0 :       break;
    1277             : 
    1278      250749 :     case t_POL:
    1279      250749 :       av = avma;
    1280      250749 :       return gerepileupto(av,poltobasis(nf,x));
    1281             : 
    1282       83668 :     case t_COL:
    1283       83668 :       if (!RgV_is_QV(x)) pari_err_TYPE("nfalgtobasis",x);
    1284       83661 :       if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
    1285       83661 :       return gcopy(x);
    1286             : 
    1287      477895 :     case t_INT:
    1288      477895 :     case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1289             :   }
    1290           0 :   pari_err_TYPE("algtobasis",x);
    1291             :   return NULL; /* LCOV_EXCL_LINE */
    1292             : }
    1293             : 
    1294             : GEN
    1295       55104 : rnfbasistoalg(GEN rnf,GEN x)
    1296             : {
    1297       55104 :   const char *f = "rnfbasistoalg";
    1298             :   long lx, i;
    1299       55104 :   pari_sp av = avma;
    1300             :   GEN z, nf, R, T;
    1301             : 
    1302       55104 :   checkrnf(rnf);
    1303       55104 :   nf = rnf_get_nf(rnf);
    1304       55104 :   T = nf_get_pol(nf);
    1305       55104 :   R = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
    1306       55104 :   switch(typ(x))
    1307             :   {
    1308         875 :     case t_COL:
    1309         875 :       z = cgetg_copy(x, &lx);
    1310        2597 :       for (i=1; i<lx; i++)
    1311             :       {
    1312        1778 :         GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
    1313        1722 :         if (typ(c) == t_POL) c = mkpolmod(c,T);
    1314        1722 :         gel(z,i) = c;
    1315             :       }
    1316         819 :       z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
    1317         735 :       return gerepileupto(av, gmodulo(z,R));
    1318             : 
    1319       34965 :     case t_POLMOD:
    1320       34965 :       x = polmod_nffix(f, rnf, x, 0);
    1321       34692 :       if (typ(x) != t_POL) break;
    1322       16046 :       retmkpolmod(RgX_copy(x), RgX_copy(R));
    1323        1582 :     case t_POL:
    1324        1582 :       if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
    1325        1337 :       if (varn(x) == varn(R))
    1326             :       {
    1327        1281 :         x = RgX_nffix(f,nf_get_pol(nf),x,0);
    1328        1281 :         return gmodulo(x, R);
    1329             :       }
    1330          56 :       pari_err_VAR(f, x,R);
    1331             :   }
    1332       36517 :   retmkpolmod(scalarpol(x, varn(R)), RgX_copy(R));
    1333             : }
    1334             : 
    1335             : GEN
    1336        2653 : matbasistoalg(GEN nf,GEN x)
    1337             : {
    1338             :   long i, j, li, lx;
    1339        2653 :   GEN z = cgetg_copy(x, &lx);
    1340             : 
    1341        2653 :   if (lx == 1) return z;
    1342        2646 :   switch(typ(x))
    1343             :   {
    1344          77 :     case t_VEC: case t_COL:
    1345         273 :       for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
    1346          77 :       return z;
    1347        2569 :     case t_MAT: break;
    1348           0 :     default: pari_err_TYPE("matbasistoalg",x);
    1349             :   }
    1350        2569 :   li = lgcols(x);
    1351        9352 :   for (j=1; j<lx; j++)
    1352             :   {
    1353        6783 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1354        6783 :     gel(z,j) = c;
    1355       30709 :     for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
    1356             :   }
    1357        2569 :   return z;
    1358             : }
    1359             : 
    1360             : GEN
    1361       31940 : matalgtobasis(GEN nf,GEN x)
    1362             : {
    1363             :   long i, j, li, lx;
    1364       31940 :   GEN z = cgetg_copy(x, &lx);
    1365             : 
    1366       31940 :   if (lx == 1) return z;
    1367       31478 :   switch(typ(x))
    1368             :   {
    1369       31471 :     case t_VEC: case t_COL:
    1370       82532 :       for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
    1371       31471 :       return z;
    1372           7 :     case t_MAT: break;
    1373           0 :     default: pari_err_TYPE("matalgtobasis",x);
    1374             :   }
    1375           7 :   li = lgcols(x);
    1376          14 :   for (j=1; j<lx; j++)
    1377             :   {
    1378           7 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1379           7 :     gel(z,j) = c;
    1380          21 :     for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
    1381             :   }
    1382           7 :   return z;
    1383             : }
    1384             : GEN
    1385       11177 : RgM_to_nfM(GEN nf,GEN x)
    1386             : {
    1387             :   long i, j, li, lx;
    1388       11177 :   GEN z = cgetg_copy(x, &lx);
    1389             : 
    1390       11177 :   if (lx == 1) return z;
    1391       11177 :   li = lgcols(x);
    1392       82810 :   for (j=1; j<lx; j++)
    1393             :   {
    1394       71633 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1395       71633 :     gel(z,j) = c;
    1396      466633 :     for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
    1397             :   }
    1398       11177 :   return z;
    1399             : }
    1400             : GEN
    1401      149392 : RgC_to_nfC(GEN nf, GEN x)
    1402      913026 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
    1403             : 
    1404             : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
    1405             : GEN
    1406      169037 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
    1407      169037 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
    1408             : GEN
    1409      169128 : polmod_nffix2(const char *f, GEN T, GEN R, GEN x, int lift)
    1410             : {
    1411      169128 :   if (RgX_equal_var(gel(x,1), R))
    1412             :   {
    1413      141036 :     x = gel(x,2);
    1414      141036 :     if (typ(x) == t_POL && varn(x) == varn(R))
    1415             :     {
    1416      106098 :       x = RgX_nffix(f, T, x, lift);
    1417      106098 :       switch(lg(x))
    1418             :       {
    1419        5831 :         case 2: return gen_0;
    1420       13604 :         case 3: return gel(x,2);
    1421             :       }
    1422       86663 :       return x;
    1423             :     }
    1424             :   }
    1425       63030 :   return Rg_nffix(f, T, x, lift);
    1426             : }
    1427             : GEN
    1428        1428 : rnfalgtobasis(GEN rnf,GEN x)
    1429             : {
    1430        1428 :   const char *f = "rnfalgtobasis";
    1431        1428 :   pari_sp av = avma;
    1432             :   GEN T, R;
    1433             : 
    1434        1428 :   checkrnf(rnf);
    1435        1428 :   R = rnf_get_pol(rnf);
    1436        1428 :   T = rnf_get_nfpol(rnf);
    1437        1428 :   switch(typ(x))
    1438             :   {
    1439          98 :     case t_COL:
    1440          98 :       if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
    1441          49 :       x = RgV_nffix(f, T, x, 0);
    1442          42 :       return gerepilecopy(av, x);
    1443             : 
    1444        1162 :     case t_POLMOD:
    1445        1162 :       x = polmod_nffix(f, rnf, x, 0);
    1446        1057 :       if (typ(x) != t_POL) break;
    1447         714 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1448         112 :     case t_POL:
    1449         112 :       if (varn(x) == varn(T))
    1450             :       {
    1451          42 :         RgX_check_QX(x,f);
    1452          28 :         if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1453          28 :         x = mkpolmod(x,T); break;
    1454             :       }
    1455          70 :       x = RgX_nffix(f, T, x, 0);
    1456          56 :       if (degpol(x) >= degpol(R)) x = RgX_rem(x, R);
    1457          56 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1458             :   }
    1459         427 :   return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
    1460             : }
    1461             : 
    1462             : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
    1463             :  * is "small" */
    1464             : GEN
    1465         259 : nfdiveuc(GEN nf, GEN a, GEN b)
    1466             : {
    1467         259 :   pari_sp av = avma;
    1468         259 :   a = nfdiv(nf,a,b);
    1469         259 :   return gerepileupto(av, ground(a));
    1470             : }
    1471             : 
    1472             : /* Given a and b in nf, gives a "small" algebraic integer r in nf
    1473             :  * of the form a-b.y */
    1474             : GEN
    1475         259 : nfmod(GEN nf, GEN a, GEN b)
    1476             : {
    1477         259 :   pari_sp av = avma;
    1478         259 :   GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
    1479         259 :   return gerepileupto(av, nfadd(nf,a,p1));
    1480             : }
    1481             : 
    1482             : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
    1483             :  * that r=a-b.y is "small". */
    1484             : GEN
    1485         259 : nfdivrem(GEN nf, GEN a, GEN b)
    1486             : {
    1487         259 :   pari_sp av = avma;
    1488         259 :   GEN p1,z, y = ground(nfdiv(nf,a,b));
    1489             : 
    1490         259 :   p1 = gneg_i(nfmul(nf,b,y));
    1491         259 :   z = cgetg(3,t_VEC);
    1492         259 :   gel(z,1) = gcopy(y);
    1493         259 :   gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
    1494             : }
    1495             : 
    1496             : /*************************************************************************/
    1497             : /**                                                                     **/
    1498             : /**                   LOGARITHMIC EMBEDDINGS                            **/
    1499             : /**                                                                     **/
    1500             : /*************************************************************************/
    1501             : 
    1502             : static int
    1503     4612157 : low_prec(GEN x)
    1504             : {
    1505     4612157 :   switch(typ(x))
    1506             :   {
    1507           0 :     case t_INT: return !signe(x);
    1508     4612157 :     case t_REAL: return !signe(x) || realprec(x) <= DEFAULTPREC;
    1509           0 :     default: return 0;
    1510             :   }
    1511             : }
    1512             : 
    1513             : static GEN
    1514       23115 : cxlog_1(GEN nf) { return zerocol(lg(nf_get_roots(nf))-1); }
    1515             : static GEN
    1516         532 : cxlog_m1(GEN nf, long prec)
    1517             : {
    1518         532 :   long i, l = lg(nf_get_roots(nf)), r1 = nf_get_r1(nf);
    1519         532 :   GEN v = cgetg(l, t_COL), p,  P;
    1520         532 :   p = mppi(prec); P = mkcomplex(gen_0, p);
    1521        1235 :   for (i = 1; i <= r1; i++) gel(v,i) = P; /* IPi*/
    1522         532 :   if (i < l) P = gmul2n(P,1);
    1523        1122 :   for (     ; i < l; i++) gel(v,i) = P; /* 2IPi */
    1524         532 :   return v;
    1525             : }
    1526             : static GEN
    1527     1715173 : ZC_cxlog(GEN nf, GEN x, long prec)
    1528             : {
    1529             :   long i, l, r1;
    1530             :   GEN v;
    1531     1715173 :   x = RgM_RgC_mul(nf_get_M(nf), Q_primpart(x));
    1532     1715173 :   l = lg(x); r1 = nf_get_r1(nf);
    1533     4330803 :   for (i = 1; i <= r1; i++)
    1534     2615630 :     if (low_prec(gel(x,i))) return NULL;
    1535     3514772 :   for (     ; i <  l;  i++)
    1536     1799599 :     if (low_prec(gnorm(gel(x,i)))) return NULL;
    1537     1715173 :   v = cgetg(l,t_COL);
    1538     4330804 :   for (i = 1; i <= r1; i++) gel(v,i) = glog(gel(x,i),prec);
    1539     3514771 :   for (     ; i <  l;  i++) gel(v,i) = gmul2n(glog(gel(x,i),prec),1);
    1540     1715172 :   return v;
    1541             : }
    1542             : static GEN
    1543      223285 : famat_cxlog(GEN nf, GEN fa, long prec)
    1544             : {
    1545      223285 :   GEN G, E, y = NULL;
    1546             :   long i, l;
    1547             : 
    1548      223285 :   if (typ(fa) != t_MAT) pari_err_TYPE("famat_cxlog",fa);
    1549      223285 :   if (lg(fa) == 1) return cxlog_1(nf);
    1550      223285 :   G = gel(fa,1);
    1551      223285 :   E = gel(fa,2); l = lg(E);
    1552     1119796 :   for (i = 1; i < l; i++)
    1553             :   {
    1554      896511 :     GEN t, e = gel(E,i), x = nf_to_scalar_or_basis(nf, gel(G,i));
    1555             :     /* multiplicative arch would be better (save logs), but exponents overflow
    1556             :      * [ could keep track of expo separately, but not worth it ] */
    1557      896511 :     switch(typ(x))
    1558             :     { /* ignore positive rationals */
    1559       16423 :       case t_FRAC: x = gel(x,1); /* fall through */
    1560      266508 :       case t_INT: if (signe(x) > 0) continue;
    1561          84 :         if (!mpodd(e)) continue;
    1562          28 :         t = cxlog_m1(nf, prec); /* we probably should not reach this line */
    1563          28 :         break;
    1564      630003 :       default: /* t_COL */
    1565      630003 :         t = ZC_cxlog(nf,x,prec); if (!t) return NULL;
    1566      630003 :         t = RgC_Rg_mul(t, e);
    1567             :     }
    1568      630031 :     y = y? RgV_add(y,t): t;
    1569             :   }
    1570      223285 :   return y ? y: cxlog_1(nf);
    1571             : }
    1572             : /* Archimedean components: [e_i Log( sigma_i(X) )], where X = primpart(x),
    1573             :  * and e_i = 1 (resp 2.) for i <= R1 (resp. > R1) */
    1574             : GEN
    1575     1309603 : nf_cxlog(GEN nf, GEN x, long prec)
    1576             : {
    1577     1309603 :   if (typ(x) == t_MAT) return famat_cxlog(nf,x,prec);
    1578     1086318 :   x = nf_to_scalar_or_basis(nf,x);
    1579     1086318 :   switch(typ(x))
    1580             :   {
    1581           0 :     case t_FRAC: x = gel(x,1); /* fall through */
    1582        1148 :     case t_INT:
    1583        1148 :       return signe(x) > 0? cxlog_1(nf): cxlog_m1(nf, prec);
    1584     1085170 :     default:
    1585     1085170 :       return ZC_cxlog(nf, x, prec);
    1586             :   }
    1587             : }
    1588             : GEN
    1589          97 : nfV_cxlog(GEN nf, GEN x, long prec)
    1590             : {
    1591             :   long i, l;
    1592          97 :   GEN v = cgetg_copy(x, &l);
    1593         167 :   for (i = 1; i < l; i++)
    1594          70 :     if (!(gel(v,i) = nf_cxlog(nf, gel(x,i), prec))) return NULL;
    1595          97 :   return v;
    1596             : }
    1597             : 
    1598             : static GEN
    1599       15239 : scalar_logembed(GEN nf, GEN u, GEN *emb)
    1600             : {
    1601             :   GEN v, logu;
    1602       15239 :   long i, s = signe(u), RU = lg(nf_get_roots(nf))-1, R1 = nf_get_r1(nf);
    1603             : 
    1604       15239 :   if (!s) pari_err_DOMAIN("nflogembed","argument","=",gen_0,u);
    1605       15239 :   v = cgetg(RU+1, t_COL); logu = logr_abs(u);
    1606       17234 :   for (i = 1; i <= R1; i++) gel(v,i) = logu;
    1607       15239 :   if (i <= RU)
    1608             :   {
    1609       14350 :     GEN logu2 = shiftr(logu,1);
    1610       55839 :     for (   ; i <= RU; i++) gel(v,i) = logu2;
    1611             :   }
    1612       15239 :   if (emb) *emb = const_col(RU, u);
    1613       15239 :   return v;
    1614             : }
    1615             : 
    1616             : static GEN
    1617        1309 : famat_logembed(GEN nf,GEN x,GEN *emb,long prec)
    1618             : {
    1619        1309 :   GEN A, M, T, a, t, g = gel(x,1), e = gel(x,2);
    1620        1309 :   long i, l = lg(e);
    1621             : 
    1622        1309 :   if (l == 1) return scalar_logembed(nf, real_1(prec), emb);
    1623        1309 :   A = NULL; T = emb? cgetg(l, t_COL): NULL;
    1624        1309 :   if (emb) *emb = M = mkmat2(T, e);
    1625       62132 :   for (i = 1; i < l; i++)
    1626             :   {
    1627       60823 :     a = nflogembed(nf, gel(g,i), &t, prec);
    1628       60823 :     if (!a) return NULL;
    1629       60823 :     a = RgC_Rg_mul(a, gel(e,i));
    1630       60823 :     A = A? RgC_add(A, a): a;
    1631       60823 :     if (emb) gel(T,i) = t;
    1632             :   }
    1633        1309 :   return A;
    1634             : }
    1635             : 
    1636             : /* Get archimedean components: [e_i log( | sigma_i(x) | )], with e_i = 1
    1637             :  * (resp 2.) for i <= R1 (resp. > R1) and set emb to the embeddings of x.
    1638             :  * Return NULL if precision problem */
    1639             : GEN
    1640       98705 : nflogembed(GEN nf, GEN x, GEN *emb, long prec)
    1641             : {
    1642             :   long i, l, r1;
    1643             :   GEN v, t;
    1644             : 
    1645       98705 :   if (typ(x) == t_MAT) return famat_logembed(nf,x,emb,prec);
    1646       97396 :   x = nf_to_scalar_or_basis(nf,x);
    1647       97396 :   if (typ(x) != t_COL) return scalar_logembed(nf, gtofp(x,prec), emb);
    1648       82157 :   x = RgM_RgC_mul(nf_get_M(nf), x);
    1649       82159 :   l = lg(x); r1 = nf_get_r1(nf); v = cgetg(l,t_COL);
    1650      109087 :   for (i = 1; i <= r1; i++)
    1651             :   {
    1652       26928 :     t = gabs(gel(x,i),prec); if (low_prec(t)) return NULL;
    1653       26928 :     gel(v,i) = glog(t,prec);
    1654             :   }
    1655      252160 :   for (   ; i < l; i++)
    1656             :   {
    1657      170002 :     t = gnorm(gel(x,i)); if (low_prec(t)) return NULL;
    1658      170002 :     gel(v,i) = glog(t,prec);
    1659             :   }
    1660       82158 :   if (emb) *emb = x;
    1661       82158 :   return v;
    1662             : }
    1663             : 
    1664             : /*************************************************************************/
    1665             : /**                                                                     **/
    1666             : /**                        REAL EMBEDDINGS                              **/
    1667             : /**                                                                     **/
    1668             : /*************************************************************************/
    1669             : static GEN
    1670      486472 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
    1671             : static GEN
    1672     1556025 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
    1673             : static GEN
    1674      608752 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
    1675             : static GEN
    1676      608753 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
    1677             : static GEN
    1678      608753 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
    1679             : 
    1680             : /* true nf, x non-zero algebraic integer; return number of positive real roots
    1681             :  * of char_x */
    1682             : static long
    1683      910167 : num_positive(GEN nf, GEN x)
    1684             : {
    1685      910167 :   GEN T = nf_get_pol(nf), B, charx;
    1686      910166 :   long dnf, vnf, N, r1 = nf_get_r1(nf);
    1687      910167 :   x = nf_to_scalar_or_alg(nf, x);
    1688      910169 :   if (typ(x) != t_POL) return (signe(x) < 0)? 0: degpol(T);
    1689             :   /* x not a scalar */
    1690      904772 :   if (r1 == 1)
    1691             :   {
    1692       31346 :     long s = signe(ZX_resultant(T, Q_primpart(x)));
    1693       31346 :     return s > 0? 1: 0;
    1694             :   }
    1695      873426 :   charx = ZXQ_charpoly(x, T, 0);
    1696      873424 :   charx = ZX_radical(charx);
    1697      873412 :   N = degpol(T) / degpol(charx);
    1698             :   /* real places are unramified ? */
    1699      873411 :   if (N == 1 || ZX_sturm(charx) * N == r1)
    1700      872815 :     return ZX_sturmpart(charx, mkvec2(gen_0,mkoo())) * N;
    1701             :   /* painful case, multiply by random square until primitive */
    1702         596 :   dnf = nf_get_degree(nf);
    1703         596 :   vnf = varn(T);
    1704         596 :   B = int2n(10);
    1705             :   for(;;)
    1706           0 :   {
    1707         596 :     GEN y = RgXQ_sqr(random_FpX(dnf, vnf, B), T);
    1708         596 :     y = RgXQ_mul(x, y, T);
    1709         596 :     charx = ZXQ_charpoly(y, T, 0);
    1710         596 :     if (ZX_is_squarefree(charx))
    1711         596 :       return ZX_sturmpart(charx, mkvec2(gen_0,mkoo()));
    1712             :   }
    1713             : }
    1714             : 
    1715             : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
    1716             :  * if x in Q. M = nf_get_M(nf) */
    1717             : static GEN
    1718        2140 : nfembed_i(GEN M, GEN x, long k)
    1719             : {
    1720        2140 :   long i, l = lg(M);
    1721        2140 :   GEN z = gel(x,1);
    1722       24380 :   for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
    1723        2140 :   return z;
    1724             : }
    1725             : GEN
    1726           0 : nfembed(GEN nf, GEN x, long k)
    1727             : {
    1728           0 :   pari_sp av = avma;
    1729           0 :   nf = checknf(nf);
    1730           0 :   x = nf_to_scalar_or_basis(nf,x);
    1731           0 :   if (typ(x) != t_COL) return gerepilecopy(av, x);
    1732           0 :   return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
    1733             : }
    1734             : 
    1735             : /* x a ZC */
    1736             : static GEN
    1737       74780 : zk_embed(GEN M, GEN x, long k)
    1738             : {
    1739       74780 :   long i, l = lg(x);
    1740       74780 :   GEN z = gel(x,1); /* times M[k,1], which is 1 */
    1741      186124 :   for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
    1742       74778 :   return z;
    1743             : }
    1744             : 
    1745             : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
    1746             : static int
    1747       24891 : oksigns(long l, GEN signs, long i, long s)
    1748             : {
    1749       24891 :   if (!signs) return s == 0;
    1750       26834 :   for (; i < l; i++)
    1751       19786 :     if (signs[i] != s) return 0;
    1752        7048 :   return 1;
    1753             : }
    1754             : 
    1755             : /* true nf, x a ZC (primitive for efficiency) which is not a scalar */
    1756             : static int
    1757       80573 : nfchecksigns_i(GEN nf, GEN x, GEN signs, GEN archp)
    1758             : {
    1759       80573 :   long i, np, npc, l = lg(archp), r1 = nf_get_r1(nf);
    1760             :   GEN sarch;
    1761             : 
    1762       80572 :   if (r1 == 0) return 1;
    1763       80179 :   np = num_positive(nf, x);
    1764       80181 :   if (np == 0)  return oksigns(l, signs, 1, 1);
    1765       71103 :   if (np == r1) return oksigns(l, signs, 1, 0);
    1766       55291 :   sarch = nfarchstar(nf, NULL, identity_perm(r1));
    1767       63850 :   for (i = 1, npc = 0; i < l; i++)
    1768             :   {
    1769       63614 :     GEN xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    1770             :     long ni, s;
    1771       63614 :     xi = Q_primpart(xi);
    1772       63614 :     ni = num_positive(nf, nfmuli(nf,x,xi));
    1773       63614 :     s = ni < np? 0: 1;
    1774       63614 :     if (s != (signs? signs[i]: 0)) return 0;
    1775       24918 :     if (!s) npc++; /* found a positive root */
    1776       24918 :     if (npc == np)
    1777             :     { /* found all positive roots */
    1778       15765 :       if (!signs) return i == l-1;
    1779        8888 :       for (i++; i < l; i++)
    1780        4236 :         if (signs[i] != 1) return 0;
    1781        4652 :       return 1;
    1782             :     }
    1783        9153 :     if (i - npc == r1 - np)
    1784             :     { /* found all negative roots */
    1785         594 :       if (!signs) return 1;
    1786         636 :       for (i++; i < l; i++)
    1787          49 :         if (signs[i]) return 0;
    1788         587 :       return 1;
    1789             :     }
    1790             :   }
    1791         236 :   return 1;
    1792             : }
    1793             : static void
    1794         985 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
    1795             : {
    1796         985 :   long i, j, l = lg(pl);
    1797         985 :   GEN signs = cgetg(l, t_VECSMALL);
    1798         985 :   GEN archp = cgetg(l, t_VECSMALL);
    1799        3080 :   for (i = j = 1; i < l; i++)
    1800             :   {
    1801        2095 :     if (!pl[i]) continue;
    1802        1578 :     archp[j] = i;
    1803        1578 :     signs[j] = (pl[i] < 0)? 1: 0;
    1804        1578 :     j++;
    1805             :   }
    1806         985 :   setlg(archp, j); *parchp = archp;
    1807         985 :   setlg(signs, j); *psigns = signs;
    1808         985 : }
    1809             : /* pl : requested signs for real embeddings, 0 = no sign constraint */
    1810             : int
    1811       15111 : nfchecksigns(GEN nf, GEN x, GEN pl)
    1812             : {
    1813       15111 :   pari_sp av = avma;
    1814             :   GEN signs, archp;
    1815       15111 :   nf = checknf(nf);
    1816       15111 :   x = nf_to_scalar_or_basis(nf,x);
    1817       15111 :   if (typ(x) != t_COL)
    1818             :   {
    1819       14126 :     long i, l = lg(pl), s = gsigne(x);
    1820       28259 :     for (i = 1; i < l; i++)
    1821       14133 :       if (pl[i] && pl[i] != s) return gc_bool(av,0);
    1822       14126 :     return gc_bool(av,1);
    1823             :   }
    1824         985 :   pl_convert(pl, &signs, &archp);
    1825         985 :   return gc_bool(av, nfchecksigns_i(nf, x, signs, archp));
    1826             : }
    1827             : 
    1828             : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
    1829             : static GEN
    1830      608754 : get_C(GEN lambda, long l, GEN signs)
    1831             : {
    1832             :   long i;
    1833             :   GEN C, mlambda;
    1834      608754 :   if (!signs) return const_vec(l-1, lambda);
    1835      579004 :   C = cgetg(l, t_COL); mlambda = gneg(lambda);
    1836     2319567 :   for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
    1837      579007 :   return C;
    1838             : }
    1839             : /* signs = NULL: totally positive at archp.
    1840             :  * Assume that a t_COL x is not a scalar */
    1841             : static GEN
    1842      722593 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
    1843             : {
    1844      722593 :   long i, l = lg(sarch_get_archp(sarch));
    1845      722596 :   GEN ex = NULL;
    1846             :   /* Is signature already correct ? */
    1847      722596 :   if (typ(x) != t_COL)
    1848             :   {
    1849      643011 :     long s = gsigne(x);
    1850      643010 :     if (!s) i = 1;
    1851      642989 :     else if (!signs)
    1852        7427 :       i = (s < 0)? 1: l;
    1853             :     else
    1854             :     {
    1855      635562 :       s = s < 0? 1: 0;
    1856     1111352 :       for (i = 1; i < l; i++)
    1857     1032608 :         if (signs[i] != s) break;
    1858             :     }
    1859      643010 :     if (i < l) ex = const_col(l-1, x);
    1860             :   }
    1861             :   else
    1862             :   { /* inefficient if x scalar, wrong if x = 0 */
    1863       79585 :     pari_sp av = avma;
    1864       79585 :     GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
    1865       79589 :     GEN xp = Q_primitive_part(x,&cex);
    1866       79588 :     if (nfchecksigns_i(nf, xp, signs, archp)) set_avma(av);
    1867             :     else
    1868             :     {
    1869       51774 :       ex = cgetg(l,t_COL);
    1870      126552 :       for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
    1871       51772 :       if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
    1872             :     }
    1873             :   }
    1874      722593 :   if (ex)
    1875             :   { /* If no, fix it */
    1876      608752 :     GEN MI = sarch_get_MI(sarch), F = sarch_get_F(sarch);
    1877      608753 :     GEN lambda = sarch_get_lambda(sarch);
    1878      608753 :     GEN t = RgC_sub(get_C(lambda, l, signs), ex);
    1879      608748 :     t = grndtoi(RgM_RgC_mul(MI,t), NULL);
    1880      608739 :     if (lg(F) != 1) t = ZM_ZC_mul(F, t);
    1881      608750 :     x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
    1882             :   }
    1883      722580 :   return x;
    1884             : }
    1885             : /* - true nf
    1886             :  * - sarch = nfarchstar(nf, F);
    1887             :  * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
    1888             :  *   (vector of signs as {0,1}-vector), NULL (totally positive at archp),
    1889             :  *   or a nonzero number field element (replaced by its signature at archp);
    1890             :  * - y is a nonzero number field element
    1891             :  * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector).
    1892             :  * Not stack-clean */
    1893             : GEN
    1894      753852 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
    1895             : {
    1896      753852 :   GEN archp = sarch_get_archp(sarch);
    1897      753853 :   if (lg(archp) == 1) return y;
    1898      720547 :   if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
    1899      720547 :   return nfsetsigns(nf, x, nf_to_scalar_or_basis(nf,y), sarch);
    1900             : }
    1901             : 
    1902             : static GEN
    1903      392037 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
    1904             : {
    1905      392037 :   GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
    1906      392044 :   lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
    1907      392037 :   if (typ(lambda) != t_REAL) lambda = gmul(lambda, uutoQ(1001,1000));
    1908      392039 :   if (lg(archp) < lg(MI))
    1909             :   {
    1910       75705 :     GEN perm = gel(indexrank(MI), 2);
    1911       75706 :     if (!F) F = matid(nf_get_degree(nf));
    1912       75706 :     MI = vecpermute(MI, perm);
    1913       75707 :     F = vecpermute(F, perm);
    1914             :   }
    1915      392043 :   if (!F) F = cgetg(1,t_MAT);
    1916      392043 :   MI = RgM_inv(MI);
    1917      392045 :   return mkvec5(DATA, archp, MI, lambda, F);
    1918             : }
    1919             : /* F nonzero integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
    1920             :  * whose sign matrix at archp is identity; archp in 'indices' format */
    1921             : GEN
    1922      567873 : nfarchstar(GEN nf, GEN F, GEN archp)
    1923             : {
    1924      567873 :   long nba = lg(archp) - 1;
    1925      567873 :   if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
    1926      389998 :   if (F && equali1(gcoeff(F,1,1))) F = NULL;
    1927      389998 :   if (F) F = idealpseudored(F, nf_get_roundG(nf));
    1928      389985 :   return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
    1929             : }
    1930             : 
    1931             : /*************************************************************************/
    1932             : /**                                                                     **/
    1933             : /**                         IDEALCHINESE                                **/
    1934             : /**                                                                     **/
    1935             : /*************************************************************************/
    1936             : static int
    1937        5305 : isprfact(GEN x)
    1938             : {
    1939             :   long i, l;
    1940             :   GEN L, E;
    1941        5305 :   if (typ(x) != t_MAT || lg(x) != 3) return 0;
    1942        5305 :   L = gel(x,1); l = lg(L);
    1943        5305 :   E = gel(x,2);
    1944       16632 :   for(i=1; i<l; i++)
    1945             :   {
    1946       11327 :     checkprid(gel(L,i));
    1947       11327 :     if (typ(gel(E,i)) != t_INT) return 0;
    1948             :   }
    1949        5305 :   return 1;
    1950             : }
    1951             : 
    1952             : /* initialize projectors mod pr[i]^e[i] for idealchinese */
    1953             : static GEN
    1954        5305 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
    1955             : {
    1956        5305 :   GEN U, E, F, FZ, L = gel(fa,1), E0 = gel(fa,2);
    1957        5305 :   long i, r = lg(L);
    1958             : 
    1959        5305 :   if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
    1960        5305 :   if (r == 1 && !dw) return cgetg(1,t_VEC);
    1961        5291 :   E = leafcopy(E0); /* do not destroy fa[2] */
    1962       16618 :   for (i = 1; i < r; i++)
    1963       11327 :     if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
    1964        5291 :   F = factorbackprime(nf, L, E);
    1965        5291 :   if (dw)
    1966             :   {
    1967         693 :     F = ZM_Z_mul(F, dw);
    1968        1596 :     for (i = 1; i < r; i++)
    1969             :     {
    1970         903 :       GEN pr = gel(L,i);
    1971         903 :       long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
    1972         903 :       if (e >= 0)
    1973         896 :         gel(E,i) = addiu(gel(E,i), v);
    1974           7 :       else if (v + e <= 0)
    1975           0 :         F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
    1976             :       else
    1977             :       {
    1978           7 :         F = idealmulpowprime(nf, F, pr, stoi(e));
    1979           7 :         gel(E,i) = stoi(v + e);
    1980             :       }
    1981             :     }
    1982             :   }
    1983        5291 :   U = cgetg(r, t_VEC);
    1984       16618 :   for (i = 1; i < r; i++)
    1985             :   {
    1986             :     GEN u;
    1987       11327 :     if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
    1988             :     else
    1989             :     {
    1990       11250 :       GEN pr = gel(L,i), e = gel(E,i), t;
    1991       11250 :       t = idealdivpowprime(nf,F, pr, e);
    1992       11250 :       u = hnfmerge_get_1(t, idealpow(nf, pr, e));
    1993       11250 :       if (!u) pari_err_COPRIME("idealchinese", t,pr);
    1994             :     }
    1995       11327 :     gel(U,i) = u;
    1996             :   }
    1997        5291 :   FZ = gcoeff(F, 1, 1);
    1998        5291 :   F = idealpseudored(F, nf_get_roundG(nf));
    1999        5291 :   return mkvec2(mkvec2(F, FZ), U);
    2000             : }
    2001             : 
    2002             : static GEN
    2003        2660 : pl_normalize(GEN nf, GEN pl)
    2004             : {
    2005        2660 :   const char *fun = "idealchinese";
    2006        2660 :   if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
    2007        2660 :   switch(typ(pl))
    2008             :   {
    2009         707 :     case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
    2010             :       /* fall through */
    2011        2660 :     case t_VECSMALL: break;
    2012           0 :     default: pari_err_TYPE(fun,pl);
    2013             :   }
    2014        2660 :   return pl;
    2015             : }
    2016             : 
    2017             : static int
    2018       11445 : is_chineseinit(GEN x)
    2019             : {
    2020             :   GEN fa, pl;
    2021             :   long l;
    2022       11445 :   if (typ(x) != t_VEC || lg(x)!=3) return 0;
    2023        9219 :   fa = gel(x,1);
    2024        9219 :   pl = gel(x,2);
    2025        9219 :   if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
    2026        5411 :   l = lg(fa);
    2027        5411 :   if (l != 1)
    2028             :   {
    2029             :     GEN z;
    2030        5369 :     if (l != 3) return 0;
    2031        5369 :     z = gel(fa, 1);
    2032        5369 :     if (typ(z) != t_VEC || lg(z) != 3 || typ(gel(z,1)) != t_MAT
    2033        5362 :                         || typ(gel(z,2)) != t_INT
    2034        5362 :                         || typ(gel(fa,2)) != t_VEC)
    2035           7 :       return 0;
    2036             :   }
    2037        5404 :   l = lg(pl);
    2038        5404 :   if (l != 1)
    2039             :   {
    2040         931 :     if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
    2041         931 :                                           || typ(gel(pl,2)) != t_VECSMALL)
    2042           0 :       return 0;
    2043             :   }
    2044        5404 :   return 1;
    2045             : }
    2046             : 
    2047             : /* nf a true 'nf' */
    2048             : static GEN
    2049        5774 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
    2050             : {
    2051        5774 :   const char *fun = "idealchineseinit";
    2052        5774 :   GEN archp = NULL, pl = NULL;
    2053        5774 :   switch(typ(fa))
    2054             :   {
    2055        2660 :     case t_VEC:
    2056        2660 :       if (is_chineseinit(fa))
    2057             :       {
    2058           0 :         if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
    2059           0 :         return fa;
    2060             :       }
    2061        2660 :       if (lg(fa) != 3) pari_err_TYPE(fun, fa);
    2062             :       /* of the form [x,s] */
    2063        2660 :       pl = pl_normalize(nf, gel(fa,2));
    2064        2660 :       fa = gel(fa,1);
    2065        2660 :       archp = vecsmall01_to_indices(pl);
    2066             :       /* keep pr_init, reset pl */
    2067        2660 :       if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
    2068             :       /* fall through */
    2069             :     case t_MAT: /* factorization? */
    2070        5305 :       if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
    2071           0 :     default: pari_err_TYPE(fun,fa);
    2072             :   }
    2073             : 
    2074        5774 :   if (!pl) pl = cgetg(1,t_VEC);
    2075             :   else
    2076             :   {
    2077        2660 :     long r = lg(archp);
    2078        2660 :     if (r == 1) pl = cgetg(1, t_VEC);
    2079             :     else
    2080             :     {
    2081        2037 :       GEN F = (lg(fa) == 1)? NULL: gmael(fa,1,1), signs = cgetg(r, t_VECSMALL);
    2082             :       long i;
    2083        5733 :       for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
    2084        2037 :       pl = setsigns_init(nf, archp, F, signs);
    2085             :     }
    2086             :   }
    2087        5774 :   return mkvec2(fa, pl);
    2088             : }
    2089             : 
    2090             : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
    2091             :  * and a vector w of elements of nf, gives b such that
    2092             :  * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
    2093             :  * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
    2094             : GEN
    2095       10709 : idealchinese(GEN nf, GEN x0, GEN w)
    2096             : {
    2097       10709 :   const char *fun = "idealchinese";
    2098       10709 :   pari_sp av = avma;
    2099       10709 :   GEN x = x0, x1, x2, s, dw, F;
    2100             : 
    2101       10709 :   nf = checknf(nf);
    2102       10709 :   if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
    2103             : 
    2104        6125 :   if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
    2105        6125 :   w = Q_remove_denom(matalgtobasis(nf,w), &dw);
    2106        6125 :   if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
    2107             :   /* x is a 'chineseinit' */
    2108        6125 :   x1 = gel(x,1); s = NULL;
    2109        6125 :   x2 = gel(x,2);
    2110        6125 :   if (lg(x1) == 1) { F = NULL; dw = NULL; }
    2111             :   else
    2112             :   {
    2113        6083 :     GEN  U = gel(x1,2), FZ;
    2114        6083 :     long i, r = lg(w);
    2115        6083 :     F = gmael(x1,1,1); FZ = gmael(x1,1,2);
    2116       20571 :     for (i=1; i<r; i++)
    2117       14488 :       if (!ZV_equal0(gel(w,i)))
    2118             :       {
    2119       10984 :         GEN t = nfmuli(nf, gel(U,i), gel(w,i));
    2120       10984 :         s = s? ZC_add(s,t): t;
    2121             :       }
    2122        6083 :     if (s)
    2123             :     {
    2124        6062 :       s = ZC_reducemodmatrix(s, F);
    2125        6062 :       if (dw && x == x0) /* input was a chineseinit */
    2126             :       {
    2127           7 :         dw = modii(dw, FZ);
    2128           7 :         s = FpC_Fp_mul(s, Fp_inv(dw, FZ), FZ);
    2129           7 :         dw = NULL;
    2130             :       }
    2131        6062 :       if (ZV_isscalar(s)) s = icopy(gel(s,1));
    2132             :     }
    2133             :   }
    2134        6125 :   if (lg(x2) != 1)
    2135             :   {
    2136        2044 :     s = nfsetsigns(nf, gel(x2,1), s? s: gen_0, x2);
    2137        2044 :     if (typ(s) == t_COL && QV_isscalar(s))
    2138             :     {
    2139         371 :       s = gel(s,1); if (!dw) s = gcopy(s);
    2140             :     }
    2141             :   }
    2142        4081 :   else if (!s) return gc_const(av, gen_0);
    2143        6076 :   return gerepileupto(av, dw? gdiv(s, dw): s);
    2144             : }
    2145             : 
    2146             : /*************************************************************************/
    2147             : /**                                                                     **/
    2148             : /**                           (Z_K/I)^*                                 **/
    2149             : /**                                                                     **/
    2150             : /*************************************************************************/
    2151             : GEN
    2152        2660 : vecsmall01_to_indices(GEN v)
    2153             : {
    2154        2660 :   long i, k, l = lg(v);
    2155        2660 :   GEN p = new_chunk(l) + l;
    2156        7567 :   for (k=1, i=l-1; i; i--)
    2157        4907 :     if (v[i]) { *--p = i; k++; }
    2158        2660 :   *--p = _evallg(k) | evaltyp(t_VECSMALL);
    2159        2660 :   set_avma((pari_sp)p); return p;
    2160             : }
    2161             : GEN
    2162     1094024 : vec01_to_indices(GEN v)
    2163             : {
    2164             :   long i, k, l;
    2165             :   GEN p;
    2166             : 
    2167     1094024 :   switch (typ(v))
    2168             :   {
    2169     1047265 :    case t_VECSMALL: return v;
    2170       46759 :    case t_VEC: break;
    2171           0 :    default: pari_err_TYPE("vec01_to_indices",v);
    2172             :   }
    2173       46759 :   l = lg(v);
    2174       46759 :   p = new_chunk(l) + l;
    2175      140584 :   for (k=1, i=l-1; i; i--)
    2176       93825 :     if (signe(gel(v,i))) { *--p = i; k++; }
    2177       46759 :   *--p = _evallg(k) | evaltyp(t_VECSMALL);
    2178       46759 :   set_avma((pari_sp)p); return p;
    2179             : }
    2180             : GEN
    2181      136891 : indices_to_vec01(GEN p, long r)
    2182             : {
    2183      136891 :   long i, l = lg(p);
    2184      136891 :   GEN v = zerovec(r);
    2185      206627 :   for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
    2186      136888 :   return v;
    2187             : }
    2188             : 
    2189             : /* return (column) vector of R1 signatures of x (0 or 1) */
    2190             : GEN
    2191     1047264 : nfsign_arch(GEN nf, GEN x, GEN arch)
    2192             : {
    2193     1047264 :   GEN sarch, V, archp = vec01_to_indices(arch);
    2194     1047265 :   long i, s, np, npc, r1, n = lg(archp)-1;
    2195             :   pari_sp av;
    2196             : 
    2197     1047265 :   if (!n) return cgetg(1,t_VECSMALL);
    2198      845122 :   if (typ(x) == t_MAT)
    2199             :   { /* factorisation */
    2200      276329 :     GEN g = gel(x,1), e = gel(x,2);
    2201      276329 :     long l = lg(g);
    2202      276329 :     V = zero_zv(n);
    2203      832041 :     for (i = 1; i < l; i++)
    2204      555709 :       if (mpodd(gel(e,i)))
    2205      435975 :         Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
    2206      276332 :     set_avma((pari_sp)V); return V;
    2207             :   }
    2208      568793 :   av = avma; V = cgetg(n+1,t_VECSMALL);
    2209      568792 :   x = nf_to_scalar_or_basis(nf, x);
    2210      568792 :   switch(typ(x))
    2211             :   {
    2212      183572 :     case t_INT:
    2213      183572 :       s = signe(x);
    2214      183572 :       if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
    2215      183572 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    2216         644 :     case t_FRAC:
    2217         644 :       s = signe(gel(x,1));
    2218         644 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    2219             :   }
    2220      384576 :   r1 = nf_get_r1(nf); x = Q_primpart(x); np = num_positive(nf, x);
    2221      384576 :   if (np == 0) { set_avma(av); return const_vecsmall(n, 1); }
    2222      338041 :   if (np == r1){ set_avma(av); return const_vecsmall(n, 0); }
    2223      253433 :   sarch = nfarchstar(nf, NULL, identity_perm(r1));
    2224      382138 :   for (i = 1, npc = 0; i <= n; i++)
    2225             :   {
    2226      381800 :     GEN xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    2227             :     long ni;
    2228      381798 :     xi = Q_primpart(xi);
    2229      381800 :     ni = num_positive(nf, nfmuli(nf,x,xi));
    2230      381800 :     V[i] = ni < np? 0: 1;
    2231      381800 :     if (!V[i]) npc++; /* found a positive root */
    2232      381800 :     if (npc == np)
    2233             :     { /* found all positive roots */
    2234      251375 :       for (i++; i <= n; i++) V[i] = 1;
    2235      136532 :       break;
    2236             :     }
    2237      245268 :     if (i - npc == r1 - np)
    2238             :     { /* found all negative roots */
    2239      181836 :       for (i++; i <= n; i++) V[i] = 0;
    2240      116565 :       break;
    2241             :     }
    2242             :   }
    2243      253435 :   set_avma((pari_sp)V); return V;
    2244             : }
    2245             : static void
    2246       36246 : chk_ind(const char *s, long i, long r1)
    2247             : {
    2248       36246 :   if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
    2249       36232 :   if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
    2250       36197 : }
    2251             : static GEN
    2252      128226 : parse_embed(GEN ind, long r, const char *f)
    2253             : {
    2254             :   long l, i;
    2255      128226 :   if (!ind) return identity_perm(r);
    2256       34111 :   switch(typ(ind))
    2257             :   {
    2258          70 :     case t_INT: ind = mkvecsmall(itos(ind)); break;
    2259          84 :     case t_VEC: case t_COL: ind = vec_to_vecsmall(ind); break;
    2260       33957 :     case t_VECSMALL: break;
    2261           0 :     default: pari_err_TYPE(f, ind);
    2262             :   }
    2263       34111 :   l = lg(ind);
    2264       70308 :   for (i = 1; i < l; i++) chk_ind(f, ind[i], r);
    2265       34062 :   return ind;
    2266             : }
    2267             : GEN
    2268      125601 : nfeltsign(GEN nf, GEN x, GEN ind0)
    2269             : {
    2270      125601 :   pari_sp av = avma;
    2271             :   long i, l;
    2272             :   GEN v, ind;
    2273      125601 :   nf = checknf(nf);
    2274      125601 :   ind = parse_embed(ind0, nf_get_r1(nf), "nfeltsign");
    2275      125580 :   l = lg(ind);
    2276      125580 :   if (is_rational_t(typ(x)))
    2277             :   { /* nfsign_arch would test this, but avoid converting t_VECSMALL -> t_VEC */
    2278             :     GEN s;
    2279       31486 :     switch(gsigne(x))
    2280             :     {
    2281       16513 :       case -1:s = gen_m1; break;
    2282       14966 :       case 1: s = gen_1; break;
    2283           7 :       default: s = gen_0; break;
    2284             :     }
    2285       31486 :     set_avma(av);
    2286       31486 :     return (ind0 && typ(ind0) == t_INT)? s: const_vec(l-1, s);
    2287             :   }
    2288       94094 :   v = nfsign_arch(nf, x, ind);
    2289       94094 :   if (ind0 && typ(ind0) == t_INT) { set_avma(av); return v[1]? gen_m1: gen_1; }
    2290       94080 :   settyp(v, t_VEC);
    2291      263928 :   for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
    2292       94080 :   return gerepileupto(av, v);
    2293             : }
    2294             : 
    2295             : /* true nf */
    2296             : GEN
    2297         728 : nfeltembed_i(GEN *pnf, GEN x, GEN ind0, long prec0)
    2298             : {
    2299             :   long i, e, l, r1, r2, prec, prec1;
    2300         728 :   GEN v, ind, cx, nf = *pnf;
    2301         728 :   nf_get_sign(nf,&r1,&r2);
    2302         728 :   x = nf_to_scalar_or_basis(nf, x);
    2303         721 :   ind = parse_embed(ind0, r1+r2, "nfeltembed");
    2304         714 :   l = lg(ind);
    2305         714 :   if (typ(x) != t_COL)
    2306             :   {
    2307         224 :     if (!(ind0 && typ(ind0) == t_INT)) x = const_vec(l-1, x);
    2308         224 :     return x;
    2309             :   }
    2310         490 :   x = Q_primitive_part(x, &cx);
    2311         490 :   prec1 = prec0; e = gexpo(x);
    2312         490 :   if (e > 8) prec1 += nbits2extraprec(e);
    2313         490 :   prec = prec1;
    2314         490 :   if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
    2315         490 :   v = cgetg(l, t_VEC);
    2316             :   for(;;)
    2317         138 :   {
    2318         628 :     GEN M = nf_get_M(nf);
    2319        2630 :     for (i = 1; i < l; i++)
    2320             :     {
    2321        2140 :       GEN t = nfembed_i(M, x, ind[i]);
    2322        2140 :       long e = gexpo(t);
    2323        2140 :       if (gequal0(t) || precision(t) < prec0
    2324        2140 :                      || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
    2325        2002 :       if (cx) t = gmul(t, cx);
    2326        2002 :       gel(v,i) = t;
    2327             :     }
    2328         628 :     if (i == l) break;
    2329         138 :     prec = precdbl(prec);
    2330         138 :     if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
    2331         138 :     *pnf = nf = nfnewprec_shallow(nf, prec);
    2332             :   }
    2333         490 :   if (ind0 && typ(ind0) == t_INT) v = gel(v,1);
    2334         490 :   return v;
    2335             : }
    2336             : GEN
    2337         728 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
    2338             : {
    2339         728 :   pari_sp av = avma; nf = checknf(nf);
    2340         728 :   return gerepilecopy(av, nfeltembed_i(&nf, x, ind0, prec0));
    2341             : }
    2342             : 
    2343             : /* number of distinct roots of sigma(f) */
    2344             : GEN
    2345        1904 : nfpolsturm(GEN nf, GEN f, GEN ind0)
    2346             : {
    2347        1904 :   pari_sp av = avma;
    2348             :   long d, l, r1, single;
    2349             :   GEN ind, u, v, vr1, T, s, t;
    2350             : 
    2351        1904 :   nf = checknf(nf); T = nf_get_pol(nf); r1 = nf_get_r1(nf);
    2352        1904 :   ind = parse_embed(ind0, r1, "nfpolsturm");
    2353        1883 :   single = ind0 && typ(ind0) == t_INT;
    2354        1883 :   l = lg(ind);
    2355             : 
    2356        1883 :   if (gequal0(f)) pari_err_ROOTS0("nfpolsturm");
    2357        1876 :   if (typ(f) == t_POL && varn(f) != varn(T))
    2358             :   {
    2359        1855 :     f = RgX_nffix("nfpolsturm", T, f,1);
    2360        1855 :     if (lg(f) == 3) f = NULL;
    2361             :   }
    2362             :   else
    2363             :   {
    2364          21 :     (void)Rg_nffix("nfpolsturm", T, f, 0);
    2365          21 :     f = NULL;
    2366             :   }
    2367        1876 :   if (!f) { set_avma(av); return single? gen_0: zerovec(l-1); }
    2368        1855 :   d = degpol(f);
    2369        1855 :   if (d == 1) { set_avma(av); return single? gen_1: const_vec(l-1,gen_1); }
    2370             : 
    2371        1785 :   vr1 = const_vecsmall(l-1, 1);
    2372        1785 :   u = Q_primpart(f); s = ZV_to_zv(nfeltsign(nf, gel(u,d+2), ind));
    2373        1785 :   v = RgX_deriv(u); t = odd(d)? leafcopy(s): zv_neg(s);
    2374             :   for(;;)
    2375         245 :   {
    2376        2030 :     GEN r = RgX_neg( Q_primpart(RgX_pseudorem(u, v)) ), sr;
    2377        2030 :     long i, dr = degpol(r);
    2378        2030 :     if (dr < 0) break;
    2379        2030 :     sr = ZV_to_zv(nfeltsign(nf, gel(r,dr+2), ind));
    2380        4865 :     for (i = 1; i < l; i++)
    2381        2835 :       if (sr[i] != s[i]) { s[i] = sr[i], vr1[i]--; }
    2382        2030 :     if (odd(dr)) sr = zv_neg(sr);
    2383        4865 :     for (i = 1; i < l; i++)
    2384        2835 :       if (sr[i] != t[i]) { t[i] = sr[i], vr1[i]++; }
    2385        2030 :     if (!dr) break;
    2386         245 :     u = v; v = r;
    2387             :   }
    2388        1785 :   if (single) return gc_stoi(av,vr1[1]);
    2389        1778 :   return gerepileupto(av, zv_to_ZV(vr1));
    2390             : }
    2391             : 
    2392             : /* True nf; return the vector of signs of x; the matrix of such if x is a vector
    2393             :  * of nf elements */
    2394             : GEN
    2395       44170 : nfsign(GEN nf, GEN x)
    2396             : {
    2397             :   long i, l;
    2398             :   GEN archp, S;
    2399             : 
    2400       44170 :   archp = identity_perm( nf_get_r1(nf) );
    2401       44170 :   if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
    2402       35945 :   l = lg(x); S = cgetg(l, t_MAT);
    2403      148111 :   for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
    2404       35945 :   return S;
    2405             : }
    2406             : 
    2407             : /* x integral elt, A integral ideal in HNF; reduce x mod A */
    2408             : static GEN
    2409     7818284 : zk_modHNF(GEN x, GEN A)
    2410     7818284 : { return (typ(x) == t_COL)?  ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
    2411             : 
    2412             : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
    2413             :    outputs an element inverse of x modulo y */
    2414             : GEN
    2415         189 : nfinvmodideal(GEN nf, GEN x, GEN y)
    2416             : {
    2417         189 :   pari_sp av = avma;
    2418         189 :   GEN a, yZ = gcoeff(y,1,1);
    2419             : 
    2420         189 :   if (equali1(yZ)) return gen_0;
    2421         189 :   x = nf_to_scalar_or_basis(nf, x);
    2422         189 :   if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
    2423             : 
    2424          79 :   a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
    2425          79 :   if (!a) pari_err_INV("nfinvmodideal", x);
    2426          79 :   return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
    2427             : }
    2428             : 
    2429             : static GEN
    2430     2688912 : nfsqrmodideal(GEN nf, GEN x, GEN id)
    2431     2688912 : { return zk_modHNF(nfsqri(nf,x), id); }
    2432             : static GEN
    2433     7292446 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
    2434     7292446 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
    2435             : /* assume x integral, k integer, A in HNF */
    2436             : GEN
    2437     5846274 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
    2438             : {
    2439     5846274 :   long s = signe(k);
    2440             :   pari_sp av;
    2441             :   GEN y;
    2442             : 
    2443     5846274 :   if (!s) return gen_1;
    2444     5846274 :   av = avma;
    2445     5846274 :   x = nf_to_scalar_or_basis(nf, x);
    2446     5846538 :   if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
    2447     2628202 :   if (s < 0) { k = negi(k); x = nfinvmodideal(nf, x,A); }
    2448     2628202 :   if (equali1(k)) return gerepileupto(av, s > 0? zk_modHNF(x, A): x);
    2449     1150505 :   for(y = NULL;;)
    2450             :   {
    2451     3839496 :     if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
    2452     3839483 :     k = shifti(k,-1); if (!signe(k)) break;
    2453     2688579 :     x = nfsqrmodideal(nf,x,A);
    2454             :   }
    2455     1150487 :   return gerepileupto(av, y);
    2456             : }
    2457             : 
    2458             : /* a * g^n mod id */
    2459             : static GEN
    2460     4695040 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
    2461             : {
    2462     4695040 :   return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
    2463             : }
    2464             : 
    2465             : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
    2466             :  * EX = multiple of exponent of (O_K/id)^* */
    2467             : GEN
    2468     2622397 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
    2469             : {
    2470     2622397 :   GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
    2471     2622397 :   long i, lx = lg(g);
    2472             : 
    2473     2622397 :   if (equali1(idZ)) return gen_1; /* id = Z_K */
    2474     2621908 :   EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
    2475     8334373 :   for (i = 1; i < lx; i++)
    2476             :   {
    2477     5712607 :     GEN h, n = centermodii(gel(e,i), EX, EXo2);
    2478     5712120 :     long sn = signe(n);
    2479     5712120 :     if (!sn) continue;
    2480             : 
    2481     4041912 :     h = nf_to_scalar_or_basis(nf, gel(g,i));
    2482     4042351 :     switch(typ(h))
    2483             :     {
    2484     2384365 :       case t_INT: break;
    2485           0 :       case t_FRAC:
    2486           0 :         h = Fp_div(gel(h,1), gel(h,2), idZ); break;
    2487     1657986 :       default:
    2488             :       {
    2489             :         GEN dh;
    2490     1657986 :         h = Q_remove_denom(h, &dh);
    2491     1658120 :         if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
    2492             :       }
    2493             :     }
    2494     4042387 :     if (sn > 0)
    2495     4040544 :       plus = nfmulpowmodideal(nf, plus, h, n, id);
    2496             :     else /* sn < 0 */
    2497        1843 :       minus = nfmulpowmodideal(nf, minus, h, negi(n), id);
    2498             :   }
    2499     2621766 :   if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
    2500     2621893 :   return plus? plus: gen_1;
    2501             : }
    2502             : 
    2503             : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
    2504             :  * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
    2505             : static GEN
    2506      237214 : zidealij(GEN x, GEN y)
    2507             : {
    2508      237214 :   GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
    2509             :   long j, N;
    2510             : 
    2511             :   /* x^(-1) y = relations between the 1 + x_i (HNF) */
    2512      237205 :   cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
    2513      237208 :   N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
    2514      574589 :   for (j=1; j<N; j++)
    2515             :   {
    2516      337407 :     GEN c = gel(G,j);
    2517      337407 :     gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
    2518      337378 :     if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
    2519             :   }
    2520      237182 :   return mkvec4(cyc, G, ZM_mul(U,xi), xp);
    2521             : }
    2522             : 
    2523             : /* lg(x) > 1, x + 1; shallow */
    2524             : static GEN
    2525      169776 : ZC_add1(GEN x)
    2526             : {
    2527      169776 :   long i, l = lg(x);
    2528      169776 :   GEN y = cgetg(l, t_COL);
    2529      396518 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2530      169778 :   gel(y,1) = addiu(gel(x,1), 1); return y;
    2531             : }
    2532             : /* lg(x) > 1, x - 1; shallow */
    2533             : static GEN
    2534       70488 : ZC_sub1(GEN x)
    2535             : {
    2536       70488 :   long i, l = lg(x);
    2537       70488 :   GEN y = cgetg(l, t_COL);
    2538      176907 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2539       70488 :   gel(y,1) = subiu(gel(x,1), 1); return y;
    2540             : }
    2541             : 
    2542             : /* x,y are t_INT or ZC */
    2543             : static GEN
    2544           0 : zkadd(GEN x, GEN y)
    2545             : {
    2546           0 :   long tx = typ(x);
    2547           0 :   if (tx == typ(y))
    2548           0 :     return tx == t_INT? addii(x,y): ZC_add(x,y);
    2549             :   else
    2550           0 :     return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
    2551             : }
    2552             : /* x a t_INT or ZC, x+1; shallow */
    2553             : static GEN
    2554      255448 : zkadd1(GEN x)
    2555             : {
    2556      255448 :   long tx = typ(x);
    2557      255448 :   return tx == t_INT? addiu(x,1): ZC_add1(x);
    2558             : }
    2559             : /* x a t_INT or ZC, x-1; shallow */
    2560             : static GEN
    2561      255486 : zksub1(GEN x)
    2562             : {
    2563      255486 :   long tx = typ(x);
    2564      255486 :   return tx == t_INT? subiu(x,1): ZC_sub1(x);
    2565             : }
    2566             : /* x,y are t_INT or ZC; x - y */
    2567             : static GEN
    2568           0 : zksub(GEN x, GEN y)
    2569             : {
    2570           0 :   long tx = typ(x), ty = typ(y);
    2571           0 :   if (tx == ty)
    2572           0 :     return tx == t_INT? subii(x,y): ZC_sub(x,y);
    2573             :   else
    2574           0 :     return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
    2575             : }
    2576             : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
    2577             : static GEN
    2578      255448 : zkmul(GEN x, GEN y)
    2579             : {
    2580      255448 :   long tx = typ(x), ty = typ(y);
    2581      255448 :   if (ty == t_INT)
    2582      184976 :     return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
    2583             :   else
    2584       70472 :     return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
    2585             : }
    2586             : 
    2587             : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
    2588             :  * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
    2589             :  * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
    2590             :  * shallow */
    2591             : GEN
    2592           0 : zkchinese(GEN zkc, GEN x, GEN y)
    2593             : {
    2594           0 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
    2595           0 :   return zk_modHNF(z, UV);
    2596             : }
    2597             : /* special case z = x mod U, = 1 mod V; shallow */
    2598             : GEN
    2599      255485 : zkchinese1(GEN zkc, GEN x)
    2600             : {
    2601      255485 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
    2602      255458 :   return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
    2603             : }
    2604             : static GEN
    2605      237464 : zkVchinese1(GEN zkc, GEN v)
    2606             : {
    2607             :   long i, ly;
    2608      237464 :   GEN y = cgetg_copy(v, &ly);
    2609      492896 :   for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
    2610      237411 :   return y;
    2611             : }
    2612             : 
    2613             : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
    2614             : GEN
    2615      237219 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
    2616             : {
    2617      237219 :   GEN v = idealaddtoone_raw(nf, A, B);
    2618             :   long e;
    2619      237203 :   if ((e = gexpo(v)) > 5)
    2620             :   {
    2621       83278 :     GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
    2622       83278 :     b= ZC_reducemodlll(b, AB);
    2623       83285 :     if (gexpo(b) < e) v = b;
    2624             :   }
    2625      237201 :   return mkvec2(zk_scalar_or_multable(nf,v), AB);
    2626             : }
    2627             : /* prepare to solve z = x (mod A), z = 1 mod (B)
    2628             :  * and then         z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
    2629             : static GEN
    2630         259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
    2631             : {
    2632         259 :   GEN zkc = zkchineseinit(nf, A, B, AB);
    2633         259 :   GEN mv = gel(zkc,1), mu;
    2634         259 :   if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
    2635          35 :   mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
    2636          35 :   return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
    2637             : }
    2638             : 
    2639             : static GEN
    2640     2156134 : apply_U(GEN L, GEN a)
    2641             : {
    2642     2156134 :   GEN e, U = gel(L,3), dU = gel(L,4);
    2643     2156134 :   if (typ(a) == t_INT)
    2644      673214 :     e = ZC_Z_mul(gel(U,1), subiu(a, 1));
    2645             :   else
    2646             :   { /* t_COL */
    2647     1482920 :     GEN t = shallowcopy(a);
    2648     1482971 :     gel(t,1) = subiu(gel(t,1), 1); /* t = a - 1 */
    2649     1482899 :     e = ZM_ZC_mul(U, t);
    2650             :   }
    2651     2156073 :   return gdiv(e, dU);
    2652             : }
    2653             : 
    2654             : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
    2655             : static GEN
    2656      169273 : principal_units(GEN nf, GEN pr, long k, GEN prk)
    2657             : {
    2658             :   GEN list, prb;
    2659      169273 :   ulong mask = quadratic_prec_mask(k);
    2660      169273 :   long a = 1;
    2661             : 
    2662      169273 :   prb = pr_hnf(nf,pr);
    2663      169278 :   list = vectrunc_init(k);
    2664      406487 :   while (mask > 1)
    2665             :   {
    2666      237216 :     GEN pra = prb;
    2667      237216 :     long b = a << 1;
    2668             : 
    2669      237216 :     if (mask & 1) b--;
    2670      237216 :     mask >>= 1;
    2671             :     /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
    2672      237216 :     prb = (b >= k)? prk: idealpows(nf,pr,b);
    2673      237216 :     vectrunc_append(list, zidealij(pra, prb));
    2674      237212 :     a = b;
    2675             :   }
    2676      169271 :   return list;
    2677             : }
    2678             : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
    2679             : static GEN
    2680     1331606 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
    2681             : {
    2682     1331606 :   GEN y = cgetg(nh+1, t_COL);
    2683     1331612 :   long j, iy, c = lg(L2)-1;
    2684     3487700 :   for (j = iy = 1; j <= c; j++)
    2685             :   {
    2686     2156133 :     GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
    2687     2155963 :     long i, nc = lg(cyc)-1;
    2688     2155963 :     int last = (j == c);
    2689     5825074 :     for (i = 1; i <= nc; i++, iy++)
    2690             :     {
    2691     3668986 :       GEN t, e = gel(E,i);
    2692     3668986 :       if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
    2693     3668979 :       t = Fp_neg(e, gel(cyc,i));
    2694     3668958 :       gel(y,iy) = negi(t);
    2695     3669110 :       if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
    2696             :     }
    2697             :   }
    2698     1331567 :   return y;
    2699             : }
    2700             : /* true nf */
    2701             : static GEN
    2702       56776 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
    2703             : {
    2704       56776 :   GEN h = cgetg(nh+1,t_MAT);
    2705       56776 :   long ih, j, c = lg(L2)-1;
    2706      181493 :   for (j = ih = 1; j <= c; j++)
    2707             :   {
    2708      124717 :     GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
    2709      124717 :     long k, lG = lg(G);
    2710      304949 :     for (k = 1; k < lG; k++,ih++)
    2711             :     { /* log(g^f) mod pr^e */
    2712      180232 :       GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
    2713      180233 :       gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
    2714      180232 :       gcoeff(h,ih,ih) = gel(F,k);
    2715             :     }
    2716             :   }
    2717       56776 :   return h;
    2718             : }
    2719             : /* true nf; k > 1; multiplicative group (1 + pr) / (1 + pr^k) */
    2720             : static GEN
    2721      169269 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
    2722             : {
    2723      169269 :   GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
    2724             : 
    2725      169273 :   L2 = principal_units(nf, pr, k, prk);
    2726      169278 :   if (k == 2)
    2727             :   {
    2728      112503 :     GEN L = gel(L2,1);
    2729      112503 :     cyc = gel(L,1);
    2730      112503 :     gen = gel(L,2);
    2731      112503 :     if (pU) *pU = matid(lg(gen)-1);
    2732             :   }
    2733             :   else
    2734             :   {
    2735       56775 :     long c = lg(L2), j;
    2736       56775 :     GEN EX, h, Ui, vg = cgetg(c, t_VEC);
    2737      181490 :     for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
    2738       56775 :     vg = shallowconcat1(vg);
    2739       56776 :     h = principal_units_relations(nf, L2, prk, lg(vg)-1);
    2740       56777 :     h = ZM_hnfall_i(h, NULL, 0);
    2741       56777 :     cyc = ZM_snf_group(h, pU, &Ui);
    2742       56777 :     c = lg(Ui); gen = cgetg(c, t_VEC); EX = cyc_get_expo(cyc);
    2743      188795 :     for (j = 1; j < c; j++)
    2744      132019 :       gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
    2745             :   }
    2746      169280 :   return mkvec4(cyc, gen, prk, L2);
    2747             : }
    2748             : GEN
    2749         182 : idealprincipalunits(GEN nf, GEN pr, long k)
    2750             : {
    2751             :   pari_sp av;
    2752             :   GEN v;
    2753         182 :   nf = checknf(nf);
    2754         182 :   if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
    2755         175 :   av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
    2756         175 :   return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
    2757             : }
    2758             : 
    2759             : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
    2760             :  * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
    2761             :  * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
    2762             :  * where
    2763             :  * cyc : type of G as abelian group (SNF)
    2764             :  * gen : generators of G, coprime to x
    2765             :  * pr^k: in HNF
    2766             :  * ff  : data for log_g in (Z_K/pr)^*
    2767             :  * Two extra components are present iff k > 1: L2, U
    2768             :  * L2  : list of data structures to compute local DL in (Z_K/pr)^*,
    2769             :  *       and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
    2770             :  * U   : base change matrices to convert a vector of local DL to DL wrt gen
    2771             :  * If MOD is not NULL, initialize G / G^MOD instead */
    2772             : static GEN
    2773      426097 : sprkinit(GEN nf, GEN pr, long k, GEN x, GEN MOD)
    2774             : {
    2775      426097 :   GEN T, p, Ld, modpr, cyc, gen, g, g0, A, prk, U, L2, ord0 = NULL;
    2776      426097 :   long f = pr_get_f(pr);
    2777             : 
    2778      426095 :   if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
    2779      426095 :   modpr = nf_to_Fq_init(nf, &pr,&T,&p);
    2780      426112 :   if (MOD)
    2781             :   {
    2782      378491 :     GEN o = subiu(powiu(p,f), 1), d = gcdii(o, MOD), fa = Z_factor(d);
    2783      378451 :     ord0 = mkvec2(o, fa); /* true order, factorization of order in G/G^MOD */
    2784      378447 :     Ld = gel(fa,1);
    2785      378447 :     if (lg(Ld) > 1 && equaliu(gel(Ld,1),2)) Ld = vecslice(Ld,2,lg(Ld)-1);
    2786             :   }
    2787             :   /* (Z_K / pr)^* */
    2788      426072 :   if (f == 1)
    2789             :   {
    2790      336886 :     g0 = g = MOD? pgener_Fp_local(p, Ld): pgener_Fp(p);
    2791      336908 :     if (!ord0) ord0 = get_arith_ZZM(subiu(p,1));
    2792             :   }
    2793             :   else
    2794             :   {
    2795       89186 :     g0 = g = MOD? gener_FpXQ_local(T, p, Ld): gener_FpXQ(T,p, &ord0);
    2796       89185 :     g = Fq_to_nf(g, modpr);
    2797       89184 :     if (typ(g) == t_POL) g = poltobasis(nf, g);
    2798             :   }
    2799      426100 :   A = gel(ord0, 1); /* Norm(pr)-1 */
    2800             :   /* If MOD != NULL, d = gcd(A, MOD): g^(A/d) has order d */
    2801      426100 :   if (k == 1)
    2802             :   {
    2803      257007 :     cyc = mkvec(A);
    2804      257005 :     gen = mkvec(g);
    2805      257001 :     prk = pr_hnf(nf,pr);
    2806      257024 :     L2 = U = NULL;
    2807             :   }
    2808             :   else
    2809             :   { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
    2810             :     GEN AB, B, u, v, w;
    2811             :     long j, l;
    2812      169093 :     w = idealprincipalunits_i(nf, pr, k, &U);
    2813             :     /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
    2814      169105 :     cyc = leafcopy(gel(w,1)); B = cyc_get_expo(cyc); AB = mulii(A,B);
    2815      169091 :     gen = leafcopy(gel(w,2));
    2816      169089 :     prk = gel(w,3);
    2817      169089 :     g = nfpowmodideal(nf, g, B, prk);
    2818      169104 :     g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
    2819      169098 :     L2 = mkvec3(A, g, gel(w,4));
    2820      169103 :     gel(cyc,1) = AB;
    2821      169103 :     gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
    2822      169085 :     u = mulii(Fp_inv(A,B), A);
    2823      169087 :     v = subui(1, u); l = lg(U);
    2824      505912 :     for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
    2825      169090 :     U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
    2826             :   }
    2827             :   /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
    2828      426121 :   if (x)
    2829             :   {
    2830      236956 :     GEN uv = zkchineseinit(nf, idealmulpowprime(nf,x,pr,utoineg(k)), prk, x);
    2831      236946 :     gen = zkVchinese1(uv, gen);
    2832             :   }
    2833      426072 :   return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
    2834             : }
    2835             : GEN
    2836     3984160 : sprk_get_cyc(GEN s) { return gel(s,1); }
    2837             : GEN
    2838     1969680 : sprk_get_expo(GEN s) { return cyc_get_expo(sprk_get_cyc(s)); }
    2839             : GEN
    2840      335933 : sprk_get_gen(GEN s) { return gel(s,2); }
    2841             : GEN
    2842     4917715 : sprk_get_prk(GEN s) { return gel(s,3); }
    2843             : GEN
    2844     2543600 : sprk_get_ff(GEN s) { return gel(s,4); }
    2845             : GEN
    2846     2604018 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
    2847             : /* L2 to 1 + pr / 1 + pr^k */
    2848             : static GEN
    2849     1213665 : sprk_get_L2(GEN s) { return gmael(s,5,3); }
    2850             : /* lift to nf of primitive root of k(pr) */
    2851             : static GEN
    2852      318219 : sprk_get_gnf(GEN s) { return gmael(s,5,2); }
    2853             : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
    2854             : void
    2855           0 : sprk_get_AgL2(GEN s, GEN *A, GEN *g, GEN *L2)
    2856           0 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
    2857             : void
    2858     1205067 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
    2859     1205067 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
    2860             : static int
    2861     2543591 : sprk_is_prime(GEN s) { return lg(s) == 5; }
    2862             : 
    2863             : GEN
    2864     1969484 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk, GEN mod)
    2865             : {
    2866     1969484 :   GEN x, expo = sprk_get_expo(sprk);
    2867     1969482 :   if (mod) expo = gcdii(expo,mod);
    2868     1969460 :   x = famat_makecoprime(nf, g, e, sprk_get_pr(sprk), sprk_get_prk(sprk), expo);
    2869     1969482 :   return log_prk(nf, x, sprk, mod);
    2870             : }
    2871             : /* famat_zlog_pr assuming (g,sprk.pr) = 1 */
    2872             : static GEN
    2873         196 : famat_zlog_pr_coprime(GEN nf, GEN g, GEN e, GEN sprk, GEN MOD)
    2874             : {
    2875         196 :   GEN x = famat_to_nf_modideal_coprime(nf, g, e, sprk_get_prk(sprk),
    2876             :                                        sprk_get_expo(sprk));
    2877         196 :   return log_prk(nf, x, sprk, MOD);
    2878             : }
    2879             : 
    2880             : /* o t_INT, O = [ord,fa] format for multiple of o (for Fq_log);
    2881             :  * return o in [ord,fa] format */
    2882             : static GEN
    2883      560267 : order_update(GEN o, GEN O)
    2884             : {
    2885      560267 :   GEN p = gmael(O,2,1), z = o, P, E;
    2886      560267 :   long i, j, l = lg(p);
    2887      560267 :   P = cgetg(l, t_COL);
    2888      560258 :   E = cgetg(l, t_COL);
    2889      617476 :   for (i = j = 1; i < l; i++)
    2890             :   {
    2891      617476 :     long v = Z_pvalrem(z, gel(p,i), &z);
    2892      617426 :     if (v)
    2893             :     {
    2894      604333 :       gel(P,j) = gel(p,i);
    2895      604333 :       gel(E,j) = utoipos(v); j++;
    2896      604357 :       if (is_pm1(z)) break;
    2897             :     }
    2898             :   }
    2899      560230 :   setlg(P, j);
    2900      560230 :   setlg(E, j); return mkvec2(o, mkmat2(P,E));
    2901             : }
    2902             : 
    2903             : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x),
    2904             :  * mod positive t_INT or NULL (meaning mod=0).
    2905             :  * return log(a) modulo mod on SNF-generators of (Z_K/pr^k)^* */
    2906             : GEN
    2907     2617543 : log_prk(GEN nf, GEN a, GEN sprk, GEN mod)
    2908             : {
    2909             :   GEN e, prk, g, U1, U2, y, ff, O, o, oN, gN,  N, T, p, modpr, pr, cyc;
    2910             : 
    2911     2617543 :   if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk, mod);
    2912     2543582 :   N = NULL;
    2913     2543582 :   ff = sprk_get_ff(sprk);
    2914     2543595 :   pr = gel(ff,1); /* modpr */
    2915     2543595 :   g = gN = gel(ff,2);
    2916     2543595 :   O = gel(ff,3); /* order of g = |Fq^*|, in [ord, fa] format */
    2917     2543595 :   o = oN = gel(O,1); /* order as a t_INT */
    2918     2543595 :   prk = sprk_get_prk(sprk);
    2919     2543606 :   modpr = nf_to_Fq_init(nf, &pr, &T, &p);
    2920     2543645 :   if (mod)
    2921             :   {
    2922     2027248 :     GEN d = gcdii(o,mod);
    2923     2026967 :     if (!equalii(o, d))
    2924             :     {
    2925      751084 :       N = diviiexact(o,d); /* > 1, coprime to p */
    2926      751029 :       a = nfpowmodideal(nf, a, N, prk);
    2927      751214 :       oN = d; /* order of g^N mod pr */
    2928             :     }
    2929             :   }
    2930     2543449 :   if (equali1(oN))
    2931      398153 :     e = gen_0;
    2932             :   else
    2933             :   {
    2934     2145369 :     if (N) { O = order_update(oN, O); gN = Fq_pow(g, N, T, p); }
    2935     2145362 :     e = Fq_log(nf_to_Fq(nf,a,modpr), gN, O, T, p);
    2936             :   }
    2937             :   /* 0 <= e < oN is correct modulo oN */
    2938     2543602 :   if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
    2939             : 
    2940      800588 :   sprk_get_U2(sprk, &U1,&U2);
    2941      800684 :   cyc = sprk_get_cyc(sprk);
    2942      800684 :   if (mod)
    2943             :   {
    2944      379326 :     cyc = ZV_snf_gcd(cyc, mod);
    2945      379314 :     if (signe(remii(mod,p))) return ZV_ZV_mod(ZC_Z_mul(U1,e), cyc);
    2946             :   }
    2947      746966 :   if (signe(e))
    2948             :   {
    2949      318219 :     GEN E = N? mulii(e, N): e;
    2950      318219 :     a = nfmulpowmodideal(nf, a, sprk_get_gnf(sprk), Fp_neg(E, o), prk);
    2951             :   }
    2952             :   /* a = 1 mod pr */
    2953      746966 :   y = log_prk1(nf, a, lg(U2)-1, sprk_get_L2(sprk), prk);
    2954      746990 :   if (N)
    2955             :   { /* from DL(a^N) to DL(a) */
    2956      135406 :     GEN E = gel(sprk_get_cyc(sprk), 1), q = powiu(p, Z_pval(E, p));
    2957      135406 :     y = ZC_Z_mul(y, Fp_inv(N, q));
    2958             :   }
    2959      746991 :   y = ZC_lincomb(gen_1, e, ZM_ZC_mul(U2,y), U1);
    2960      746990 :   return ZV_ZV_mod(y, cyc);
    2961             : }
    2962             : /* true nf */
    2963             : GEN
    2964       90237 : log_prk_init(GEN nf, GEN pr, long k, GEN MOD)
    2965       90237 : { return sprkinit(nf,pr,k,NULL,MOD);}
    2966             : GEN
    2967         497 : veclog_prk(GEN nf, GEN v, GEN sprk)
    2968             : {
    2969         497 :   long l = lg(v), i;
    2970         497 :   GEN w = cgetg(l, t_MAT);
    2971        1232 :   for (i = 1; i < l; i++) gel(w,i) = log_prk(nf, gel(v,i), sprk, NULL);
    2972         497 :   return w;
    2973             : }
    2974             : 
    2975             : static GEN
    2976     1374185 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
    2977             : {
    2978     1374185 :   long i, l0, l = lg(S->U);
    2979     1374185 :   GEN g = gel(fa,1), e = gel(fa,2), y = cgetg(l, t_COL);
    2980     1374186 :   l0 = lg(S->sprk); /* = l (trivial arch. part), or l-1 */
    2981     2852230 :   for (i=1; i < l0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i), S->mod);
    2982     1374185 :   if (l0 != l)
    2983             :   {
    2984      190902 :     if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
    2985      190902 :     gel(y,l0) = Flc_to_ZC(sgn);
    2986             :   }
    2987     1374185 :   return y;
    2988             : }
    2989             : 
    2990             : /* assume that cyclic factors are normalized, in particular != [1] */
    2991             : static GEN
    2992      257545 : split_U(GEN U, GEN Sprk)
    2993             : {
    2994      257545 :   long t = 0, k, n, l = lg(Sprk);
    2995      257545 :   GEN vU = cgetg(l+1, t_VEC);
    2996      592703 :   for (k = 1; k < l; k++)
    2997             :   {
    2998      335155 :     n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
    2999      335156 :     gel(vU,k) = vecslice(U, t+1, t+n);
    3000      335159 :     t += n;
    3001             :   }
    3002             :   /* t+1 .. lg(U)-1 */
    3003      257548 :   n = lg(U) - t - 1; /* can be 0 */
    3004      257548 :   if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
    3005      257559 :   return vU;
    3006             : }
    3007             : 
    3008             : static void
    3009     1990759 : init_zlog_mod(zlog_S *S, GEN bid, GEN mod)
    3010             : {
    3011     1990759 :   GEN fa2 = bid_get_fact2(bid), MOD = bid_get_MOD(bid);
    3012     1990748 :   S->U = bid_get_U(bid);
    3013     1990742 :   S->hU = lg(bid_get_cyc(bid))-1;
    3014     1990734 :   S->archp = bid_get_archp(bid);
    3015     1990736 :   S->sprk = bid_get_sprk(bid);
    3016     1990737 :   S->bid = bid;
    3017     1990737 :   if (MOD) mod = mod? gcdii(mod, MOD): MOD;
    3018     1990637 :   S->mod = mod;
    3019     1990637 :   S->P = gel(fa2,1);
    3020     1990637 :   S->k = gel(fa2,2);
    3021     1990637 :   S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
    3022     1990662 : }
    3023             : void
    3024      380272 : init_zlog(zlog_S *S, GEN bid)
    3025             : {
    3026      380272 :   return init_zlog_mod(S, bid, NULL);
    3027             : }
    3028             : 
    3029             : /* a a t_FRAC/t_INT, reduce mod bid */
    3030             : static GEN
    3031          14 : Q_mod_bid(GEN bid, GEN a)
    3032             : {
    3033          14 :   GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
    3034          14 :   GEN b = Rg_to_Fp(a, xZ);
    3035          14 :   if (gsigne(a) < 0) b = subii(b, xZ);
    3036          14 :   return signe(b)? b: xZ;
    3037             : }
    3038             : /* Return decomposition of a on the CRT generators blocks attached to the
    3039             :  * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
    3040             : static GEN
    3041      381548 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
    3042             : {
    3043             :   long k, l;
    3044             :   GEN y;
    3045      381548 :   a = nf_to_scalar_or_basis(nf, a);
    3046      381543 :   switch(typ(a))
    3047             :   {
    3048      162588 :     case t_INT: break;
    3049          14 :     case t_FRAC: a = Q_mod_bid(S->bid, a); break;
    3050      218941 :     default: /* case t_COL: */
    3051             :     {
    3052             :       GEN den;
    3053      218941 :       check_nfelt(a, &den);
    3054      218956 :       if (den)
    3055             :       {
    3056         105 :         a = Q_muli_to_int(a, den);
    3057         105 :         a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
    3058         105 :         return famat_zlog(nf, a, sgn, S);
    3059             :       }
    3060             :     }
    3061             :   }
    3062      381448 :   if (sgn)
    3063      374546 :     sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
    3064             :   else
    3065        6902 :     sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
    3066      381452 :   l = lg(S->sprk);
    3067      381452 :   y = cgetg(sgn? l+1: l, t_COL);
    3068      922807 :   for (k = 1; k < l; k++)
    3069             :   {
    3070      541422 :     GEN sprk = gel(S->sprk,k);
    3071      541422 :     gel(y,k) = log_prk(nf, a, sprk, S->mod);
    3072             :   }
    3073      381385 :   if (sgn) gel(y,l) = Flc_to_ZC(sgn);
    3074      381397 :   return y;
    3075             : }
    3076             : 
    3077             : /* true nf */
    3078             : GEN
    3079       43813 : pr_basis_perm(GEN nf, GEN pr)
    3080             : {
    3081       43813 :   long f = pr_get_f(pr);
    3082             :   GEN perm;
    3083       43813 :   if (f == nf_get_degree(nf)) return identity_perm(f);
    3084       38164 :   perm = cgetg(f+1, t_VECSMALL);
    3085       38164 :   perm[1] = 1;
    3086       38164 :   if (f > 1)
    3087             :   {
    3088        2912 :     GEN H = pr_hnf(nf,pr);
    3089             :     long i, k;
    3090       10808 :     for (i = k = 2; k <= f; i++)
    3091        7896 :       if (!equali1(gcoeff(H,i,i))) perm[k++] = i;
    3092             :   }
    3093       38164 :   return perm;
    3094             : }
    3095             : 
    3096             : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
    3097             : static GEN
    3098     1755688 : ZMV_ZCV_mul(GEN U, GEN y)
    3099             : {
    3100     1755688 :   long i, l = lg(U);
    3101     1755688 :   GEN z = NULL;
    3102     1755688 :   if (l == 1) return cgetg(1,t_COL);
    3103     4140106 :   for (i = 1; i < l; i++)
    3104             :   {
    3105     2384501 :     GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
    3106     2384465 :     z = z? ZC_add(z, u): u;
    3107             :   }
    3108     1755605 :   return z;
    3109             : }
    3110             : /* A * (U[1], ..., U[d] */
    3111             : static GEN
    3112         518 : ZM_ZMV_mul(GEN A, GEN U)
    3113             : {
    3114             :   long i, l;
    3115         518 :   GEN V = cgetg_copy(U,&l);
    3116        1057 :   for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
    3117         518 :   return V;
    3118             : }
    3119             : 
    3120             : /* a = 1 mod pr, sprk mod pr^e, e >= 1 */
    3121             : static GEN
    3122      404408 : sprk_log_prk1_2(GEN nf, GEN a, GEN sprk)
    3123             : {
    3124      404408 :   GEN U1, U2, y, L2 = sprk_get_L2(sprk);
    3125      404408 :   sprk_get_U2(sprk, &U1,&U2);
    3126      404407 :   y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, sprk_get_prk(sprk)));
    3127      404394 :   return ZV_ZV_mod(y, sprk_get_cyc(sprk));
    3128             : }
    3129             : /* true nf; assume e >= 2 */
    3130             : GEN
    3131      105867 : sprk_log_gen_pr2(GEN nf, GEN sprk, long e)
    3132             : {
    3133      105867 :   GEN M, G, pr = sprk_get_pr(sprk);
    3134             :   long i, l;
    3135      105866 :   if (e == 2)
    3136             :   {
    3137       62307 :     GEN L2 = sprk_get_L2(sprk), L = gel(L2,1);
    3138       62307 :     G = gel(L,2); l = lg(G);
    3139             :   }
    3140             :   else
    3141             :   {
    3142       43559 :     GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
    3143       43560 :     l = lg(perm);
    3144       43560 :     G = cgetg(l, t_VEC);
    3145       43561 :     if (typ(PI) == t_INT)
    3146             :     { /* zk_ei_mul doesn't allow t_INT */
    3147        5642 :       long N = nf_get_degree(nf);
    3148        5642 :       gel(G,1) = addiu(PI,1);
    3149        8645 :       for (i = 2; i < l; i++)
    3150             :       {
    3151        3003 :         GEN z = col_ei(N, 1);
    3152        3003 :         gel(G,i) = z; gel(z, perm[i]) = PI;
    3153             :       }
    3154             :     }
    3155             :     else
    3156             :     {
    3157       37919 :       gel(G,1) = nfadd(nf, gen_1, PI);
    3158       44702 :       for (i = 2; i < l; i++)
    3159        6783 :         gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
    3160             :     }
    3161             :   }
    3162      105868 :   M = cgetg(l, t_MAT);
    3163      234396 :   for (i = 1; i < l; i++) gel(M,i) = sprk_log_prk1_2(nf, gel(G,i), sprk);
    3164      105844 :   return M;
    3165             : }
    3166             : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
    3167             :  * defined implicitly via CRT. 'ind' is the index of pr in modulus
    3168             :  * factorization; true nf */
    3169             : GEN
    3170      413968 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
    3171             : {
    3172      413968 :   GEN Uind = gel(S->U, ind);
    3173      413968 :   if (e == 1) retmkmat( gel(Uind,1) );
    3174      103164 :   return ZM_mul(Uind, sprk_log_gen_pr2(nf, gel(S->sprk,ind), e));
    3175             : }
    3176             : /* true nf */
    3177             : GEN
    3178        2037 : sprk_log_gen_pr(GEN nf, GEN sprk, long e)
    3179             : {
    3180        2037 :   if (e == 1)
    3181             :   {
    3182           0 :     long n = lg(sprk_get_cyc(sprk))-1;
    3183           0 :     retmkmat(col_ei(n, 1));
    3184             :   }
    3185        2037 :   return sprk_log_gen_pr2(nf, sprk, e);
    3186             : }
    3187             : /* a = 1 mod pr */
    3188             : GEN
    3189      275855 : sprk_log_prk1(GEN nf, GEN a, GEN sprk)
    3190             : {
    3191      275855 :   if (lg(sprk) == 5) return mkcol(gen_0); /* mod pr */
    3192      275855 :   return sprk_log_prk1_2(nf, a, sprk);
    3193             : }
    3194             : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
    3195             :  * v = vector of r1 real places */
    3196             : GEN
    3197       86265 : log_gen_arch(zlog_S *S, long index) { return gel(veclast(S->U), index); }
    3198             : 
    3199             : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
    3200             : static GEN
    3201      258573 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
    3202             : {
    3203      258573 :   GEN G, h = ZV_prod(cyc);
    3204             :   long c;
    3205      258591 :   if (!U) return mkvec2(h,cyc);
    3206      258234 :   c = lg(U);
    3207      258234 :   G = cgetg(c,t_VEC);
    3208      258243 :   if (c > 1)
    3209             :   {
    3210      228146 :     GEN U0, Uoo, EX = cyc_get_expo(cyc); /* exponent of bid */
    3211      228145 :     long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
    3212      228155 :     if (!nba) { U0 = U; Uoo = NULL; }
    3213       80419 :     else if (nba == hU) { U0 = NULL; Uoo = U; }
    3214             :     else
    3215             :     {
    3216       71277 :       U0 = rowslice(U, 1, hU-nba);
    3217       71280 :       Uoo = rowslice(U, hU-nba+1, hU);
    3218             :     }
    3219      695681 :     for (i = 1; i < c; i++)
    3220             :     {
    3221      467545 :       GEN t = gen_1;
    3222      467545 :       if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
    3223      467533 :       if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
    3224      467524 :       gel(G,i) = t;
    3225             :     }
    3226             :   }
    3227      258233 :   return mkvec3(h, cyc, G);
    3228             : }
    3229             : 
    3230             : /* remove prime ideals of norm 2 with exponent 1 from factorization */
    3231             : static GEN
    3232      258924 : famat_strip2(GEN fa)
    3233             : {
    3234      258924 :   GEN P = gel(fa,1), E = gel(fa,2), Q, F;
    3235      258924 :   long l = lg(P), i, j;
    3236      258924 :   Q = cgetg(l, t_COL);
    3237      258919 :   F = cgetg(l, t_COL);
    3238      634099 :   for (i = j = 1; i < l; i++)
    3239             :   {
    3240      375176 :     GEN pr = gel(P,i), e = gel(E,i);
    3241      375176 :     if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
    3242             :     {
    3243      336545 :       gel(Q,j) = pr;
    3244      336545 :       gel(F,j) = e; j++;
    3245             :     }
    3246             :   }
    3247      258923 :   setlg(Q,j);
    3248      258922 :   setlg(F,j); return mkmat2(Q,F);
    3249             : }
    3250             : static int
    3251      134091 : checkarchp(GEN v, long r1)
    3252             : {
    3253      134091 :   long i, l = lg(v);
    3254      134091 :   pari_sp av = avma;
    3255             :   GEN p;
    3256      134091 :   if (l == 1) return 1;
    3257       47156 :   if (l == 2) return v[1] > 0 && v[1] <= r1;
    3258       22019 :   p = zero_zv(r1);
    3259       66150 :   for (i = 1; i < l; i++)
    3260             :   {
    3261       44163 :     long j = v[i];
    3262       44163 :     if (j <= 0 || j > r1 || p[j]) return gc_long(av, 0);
    3263       44128 :     p[j] = 1;
    3264             :   }
    3265       21987 :   return gc_long(av, 1);
    3266             : }
    3267             : 
    3268             : /* True nf. Put ideal to form [[ideal,arch]] and set fa and fa2 to its
    3269             :  * factorization, archp to the indices of arch places */
    3270             : GEN
    3271      258920 : check_mod_factored(GEN nf, GEN ideal, GEN *fa_, GEN *fa2_, GEN *archp_, GEN MOD)
    3272             : {
    3273             :   GEN arch, x, fa, fa2, archp;
    3274             :   long R1;
    3275             : 
    3276      258920 :   R1 = nf_get_r1(nf);
    3277      258919 :   if (typ(ideal) == t_VEC && lg(ideal) == 3)
    3278             :   {
    3279      178715 :     arch = gel(ideal,2);
    3280      178715 :     ideal= gel(ideal,1);
    3281      178715 :     switch(typ(arch))
    3282             :     {
    3283       44624 :       case t_VEC:
    3284       44624 :         if (lg(arch) != R1+1)
    3285           7 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3286       44617 :         archp = vec01_to_indices(arch);
    3287       44617 :         break;
    3288      134091 :       case t_VECSMALL:
    3289      134091 :         if (!checkarchp(arch, R1))
    3290          35 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3291      134054 :         archp = arch;
    3292      134054 :         arch = indices_to_vec01(archp, R1);
    3293      134053 :         break;
    3294           0 :       default:
    3295           0 :         pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    3296             :         return NULL;/*LCOV_EXCL_LINE*/
    3297             :     }
    3298             :   }
    3299             :   else
    3300             :   {
    3301       80204 :     arch = zerovec(R1);
    3302       80206 :     archp = cgetg(1, t_VECSMALL);
    3303             :   }
    3304      258875 :   if (MOD)
    3305             :   {
    3306      214271 :     if (typ(MOD) != t_INT) pari_err_TYPE("bnrinit [incorrect cycmod]", MOD);
    3307      214271 :     if (mpodd(MOD) && lg(archp) != 1)
    3308         231 :       MOD = shifti(MOD, 1); /* ensure elements of G^MOD are >> 0 */
    3309             :   }
    3310      258874 :   if (is_nf_factor(ideal))
    3311             :   {
    3312       40362 :     fa = ideal;
    3313       40362 :     x = factorbackprime(nf, gel(fa,1), gel(fa,2));
    3314             :   }
    3315             :   else
    3316             :   {
    3317      218510 :     fa = idealfactor(nf, ideal);
    3318      218524 :     x = ideal;
    3319             :   }
    3320      258886 :   if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
    3321      258868 :   if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
    3322      258868 :   if (typ(gcoeff(x,1,1)) != t_INT)
    3323           7 :     pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
    3324             : 
    3325      258861 :   fa2 = famat_strip2(fa);
    3326      258851 :   if (fa_ != NULL) *fa_ = fa;
    3327      258851 :   if (fa2_ != NULL) *fa2_ = fa2;
    3328      258851 :   if (fa2_ != NULL) *archp_ = archp;
    3329      258851 :   return mkvec2(x, arch);
    3330             : }
    3331             : 
    3332             : /* True nf. Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
    3333             :    flag may include nf_GEN | nf_INIT */
    3334             : static GEN
    3335      258281 : Idealstarmod_i(GEN nf, GEN ideal, long flag, GEN MOD)
    3336             : {
    3337             :   long i, nbp;
    3338      258281 :   GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x_arch, x, arch, archp, E, P, sarch, gen;
    3339             : 
    3340      258281 :   x_arch = check_mod_factored(nf, ideal, &fa, &fa2, &archp, MOD);
    3341      258219 :   x = gel(x_arch, 1);
    3342      258219 :   arch = gel(x_arch, 2);
    3343             : 
    3344      258219 :   sarch = nfarchstar(nf, x, archp);
    3345      258224 :   P = gel(fa2,1);
    3346      258224 :   E = gel(fa2,2);
    3347      258224 :   nbp = lg(P)-1;
    3348      258224 :   sprk = cgetg(nbp+1,t_VEC);
    3349      258231 :   if (nbp)
    3350             :   {
    3351      218949 :     GEN t = (lg(gel(fa,1))==2)? NULL: x; /* beware fa != fa2 */
    3352      218949 :     cyc = cgetg(nbp+2,t_VEC);
    3353      218937 :     gen = cgetg(nbp+1,t_VEC);
    3354      554822 :     for (i = 1; i <= nbp; i++)
    3355             :     {
    3356      335855 :       GEN L = sprkinit(nf, gel(P,i), itou(gel(E,i)), t, MOD);
    3357      335870 :       gel(sprk,i) = L;
    3358      335870 :       gel(cyc,i) = sprk_get_cyc(L);
    3359             :       /* true gens are congruent to those mod x AND positive at archp */
    3360      335870 :       gel(gen,i) = sprk_get_gen(L);
    3361             :     }
    3362      218967 :     gel(cyc,i) = sarch_get_cyc(sarch);
    3363      218966 :     cyc = shallowconcat1(cyc);
    3364      218969 :     gen = shallowconcat1(gen);
    3365      218973 :     cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
    3366             :   }
    3367             :   else
    3368             :   {
    3369       39282 :     cyc = sarch_get_cyc(sarch);
    3370       39282 :     gen = cgetg(1,t_VEC);
    3371       39282 :     U = matid(lg(cyc)-1);
    3372       39282 :     if (flag & nf_GEN) u1 = U;
    3373             :   }
    3374      258230 :   if (MOD) cyc = ZV_snf_gcd(cyc, MOD);
    3375      258211 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3376      258238 :   if (!(flag & nf_INIT)) return y;
    3377      257440 :   U = split_U(U, sprk);
    3378      514892 :   return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2),
    3379      257445 :                 MOD? mkvec3(sprk, sarch, MOD): mkvec2(sprk, sarch),
    3380             :                 U);
    3381             : }
    3382             : 
    3383             : static long
    3384          63 : idealHNF_norm_pval(GEN x, GEN p)
    3385             : {
    3386          63 :   long i, v = 0, l = lg(x);
    3387         175 :   for (i = 1; i < l; i++) v += Z_pval(gcoeff(x,i,i), p);
    3388          63 :   return v;
    3389             : }
    3390             : static long
    3391          63 : sprk_get_k(GEN sprk)
    3392             : {
    3393             :   GEN pr, prk;
    3394          63 :   if (sprk_is_prime(sprk)) return 1;
    3395          63 :   pr = sprk_get_pr(sprk);
    3396          63 :   prk = sprk_get_prk(sprk);
    3397          63 :   return idealHNF_norm_pval(prk, pr_get_p(pr)) / pr_get_f(pr);
    3398             : }
    3399             : /* true nf, L a sprk */
    3400             : GEN
    3401          63 : sprk_to_bid(GEN nf, GEN L, long flag)
    3402             : {
    3403          63 :   GEN y, cyc, U, u1 = NULL, fa, fa2, arch, sarch, gen, sprk;
    3404             : 
    3405          63 :   arch = zerovec(nf_get_r1(nf));
    3406          63 :   fa = to_famat_shallow(sprk_get_pr(L), utoipos(sprk_get_k(L)));
    3407          63 :   sarch = nfarchstar(nf, NULL, cgetg(1, t_VECSMALL));
    3408          63 :   fa2 = famat_strip2(fa);
    3409          63 :   sprk = mkvec(L);
    3410          63 :   cyc = shallowconcat(sprk_get_cyc(L), sarch_get_cyc(sarch));
    3411          63 :   gen = sprk_get_gen(L);
    3412          63 :   cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
    3413          63 :   y = bid_grp(nf, u1, cyc, gen, NULL, sarch);
    3414          63 :   if (!(flag & nf_INIT)) return y;
    3415          63 :   return mkvec5(mkvec2(sprk_get_prk(L), arch), y, mkvec2(fa,fa2),
    3416             :                 mkvec2(sprk, sarch), split_U(U, sprk));
    3417             : }
    3418             : GEN
    3419      258011 : Idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
    3420             : {
    3421      258011 :   pari_sp av = avma;
    3422      258011 :   nf = nf? checknf(nf): nfinit(pol_x(0), DEFAULTPREC);
    3423      258008 :   return gerepilecopy(av, Idealstarmod_i(nf, ideal, flag, MOD));
    3424             : }
    3425             : GEN
    3426         938 : Idealstar(GEN nf, GEN ideal, long flag) { return Idealstarmod(nf, ideal, flag, NULL); }
    3427             : GEN
    3428         273 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
    3429             : {
    3430         273 :   pari_sp av = avma;
    3431         273 :   GEN z = Idealstarmod_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag, NULL);
    3432         273 :   return gerepilecopy(av, z);
    3433             : }
    3434             : 
    3435             : /* FIXME: obsolete */
    3436             : GEN
    3437           0 : zidealstarinitgen(GEN nf, GEN ideal)
    3438           0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
    3439             : GEN
    3440           0 : zidealstarinit(GEN nf, GEN ideal)
    3441           0 : { return Idealstar(nf,ideal, nf_INIT); }
    3442             : GEN
    3443           0 : zidealstar(GEN nf, GEN ideal)
    3444           0 : { return Idealstar(nf,ideal, nf_GEN); }
    3445             : 
    3446             : GEN
    3447         112 : idealstarmod(GEN nf, GEN ideal, long flag, GEN MOD)
    3448             : {
    3449         112 :   switch(flag)
    3450             :   {
    3451           0 :     case 0: return Idealstarmod(nf,ideal, nf_GEN, MOD);
    3452          98 :     case 1: return Idealstarmod(nf,ideal, nf_INIT, MOD);
    3453          14 :     case 2: return Idealstarmod(nf,ideal, nf_INIT|nf_GEN, MOD);
    3454           0 :     default: pari_err_FLAG("idealstar");
    3455             :   }
    3456             :   return NULL; /* LCOV_EXCL_LINE */
    3457             : }
    3458             : GEN
    3459           0 : idealstar0(GEN nf, GEN ideal,long flag) { return idealstarmod(nf, ideal, flag, NULL); }
    3460             : 
    3461             : void
    3462      218953 : check_nfelt(GEN x, GEN *den)
    3463             : {
    3464      218953 :   long l = lg(x), i;
    3465      218953 :   GEN t, d = NULL;
    3466      218953 :   if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
    3467      809129 :   for (i=1; i<l; i++)
    3468             :   {
    3469      590174 :     t = gel(x,i);
    3470      590174 :     switch (typ(t))
    3471             :     {
    3472      589943 :       case t_INT: break;
    3473         231 :       case t_FRAC:
    3474         231 :         if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
    3475         231 :         break;
    3476           0 :       default: pari_err_TYPE("check_nfelt", x);
    3477             :     }
    3478             :   }
    3479      218955 :   *den = d;
    3480      218955 : }
    3481             : 
    3482             : GEN
    3483     1953073 : ZV_snf_gcd(GEN x, GEN mod)
    3484     4358611 : { pari_APPLY_same(gcdii(gel(x,i), mod)); }
    3485             : 
    3486             : /* assume a true bnf and bid */
    3487             : GEN
    3488      227126 : ideallog_units0(GEN bnf, GEN bid, GEN MOD)
    3489             : {
    3490      227126 :   GEN nf = bnf_get_nf(bnf), D, y, C, cyc;
    3491      227124 :   long j, lU = lg(bnf_get_logfu(bnf)); /* r1+r2 */
    3492             :   zlog_S S;
    3493      227124 :   init_zlog_mod(&S, bid, MOD);
    3494      227104 :   if (!S.hU) return zeromat(0,lU);
    3495      227104 :   cyc = bid_get_cyc(bid);
    3496      227102 :   D = nfsign_fu(bnf, bid_get_archp(bid));
    3497      227113 :   y = cgetg(lU, t_MAT);
    3498      227113 :   if ((C = bnf_build_cheapfu(bnf)))
    3499      374496 :   { for (j = 1; j < lU; j++) gel(y,j) = zlog(nf, gel(C,j), gel(D,j), &S); }
    3500             :   else
    3501             :   {
    3502          49 :     long i, l = lg(S.U), l0 = lg(S.sprk);
    3503             :     GEN X, U;
    3504          49 :     if (!(C = bnf_compactfu_mat(bnf))) bnf_build_units(bnf); /* error */
    3505          49 :     X = gel(C,1); U = gel(C,2);
    3506         147 :     for (j = 1; j < lU; j++) gel(y,j) = cgetg(l, t_COL);
    3507         126 :     for (i = 1; i < l0; i++)
    3508             :     {
    3509          77 :       GEN sprk = gel(S.sprk, i);
    3510          77 :       GEN Xi = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3511         231 :       for (j = 1; j < lU; j++)
    3512         154 :         gcoeff(y,i,j) = famat_zlog_pr_coprime(nf, Xi, gel(U,j), sprk, MOD);
    3513             :     }
    3514          49 :     if (l0 != l)
    3515          56 :       for (j = 1; j < lU; j++) gcoeff(y,l0,j) = Flc_to_ZC(gel(D,j));
    3516             :   }
    3517      227112 :   y = vec_prepend(y, zlog(nf, bnf_get_tuU(bnf), nfsign_tu(bnf, S.archp), &S));
    3518      601721 :   for (j = 1; j <= lU; j++)
    3519      374616 :     gel(y,j) = ZV_ZV_mod(ZMV_ZCV_mul(S.U, gel(y,j)), cyc);
    3520      227105 :   return y;
    3521             : }
    3522             : GEN
    3523          84 : ideallog_units(GEN bnf, GEN bid)
    3524          84 : { return ideallog_units0(bnf, bid, NULL); }
    3525             : GEN
    3526         518 : log_prk_units(GEN nf, GEN D, GEN sprk)
    3527             : {
    3528         518 :   GEN L, Ltu = log_prk(nf, gel(D,1), sprk, NULL);
    3529         518 :   D = gel(D,2);
    3530         518 :   if (lg(D) != 3 || typ(gel(D,2)) != t_MAT) L = veclog_prk(nf, D, sprk);
    3531             :   else
    3532             :   {
    3533          21 :     GEN X = gel(D,1), U = gel(D,2);
    3534          21 :     long j, lU = lg(U);
    3535          21 :     X = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3536          21 :     L = cgetg(lU, t_MAT);
    3537          63 :     for (j = 1; j < lU; j++)
    3538          42 :       gel(L,j) = famat_zlog_pr_coprime(nf, X, gel(U,j), sprk, NULL);
    3539             :   }
    3540         518 :   return vec_prepend(L, Ltu);
    3541             : }
    3542             : 
    3543             : static GEN
    3544     1383381 : ideallog_i(GEN nf, GEN x, zlog_S *S)
    3545             : {
    3546     1383381 :   pari_sp av = avma;
    3547             :   GEN y;
    3548     1383381 :   if (!S->hU) return cgetg(1, t_COL);
    3549     1381085 :   if (typ(x) == t_MAT)
    3550     1374078 :     y = famat_zlog(nf, x, NULL, S);
    3551             :   else
    3552        7007 :     y = zlog(nf, x, NULL, S);
    3553     1381082 :   y = ZMV_ZCV_mul(S->U, y);
    3554     1381078 :   return gerepileupto(av, ZV_ZV_mod(y, bid_get_cyc(S->bid)));
    3555             : }
    3556             : GEN
    3557     1390062 : ideallogmod(GEN nf, GEN x, GEN bid, GEN mod)
    3558             : {
    3559             :   zlog_S S;
    3560     1390062 :   if (!nf)
    3561             :   {
    3562        6671 :     if (mod) pari_err_IMPL("Zideallogmod");
    3563        6671 :     return Zideallog(bid, x);
    3564             :   }
    3565     1383391 :   checkbid(bid); init_zlog_mod(&S, bid, mod);
    3566     1383381 :   return ideallog_i(checknf(nf), x, &S);
    3567             : }
    3568             : GEN
    3569       13769 : ideallog(GEN nf, GEN x, GEN bid) { return ideallogmod(nf, x, bid, NULL); }
    3570             : 
    3571             : /*************************************************************************/
    3572             : /**                                                                     **/
    3573             : /**               JOIN BID STRUCTURES, IDEAL LISTS                      **/
    3574             : /**                                                                     **/
    3575             : /*************************************************************************/
    3576             : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
    3577             :  * Output: bid for m1 m2 */
    3578             : static GEN
    3579         469 : join_bid(GEN nf, GEN bid1, GEN bid2)
    3580             : {
    3581         469 :   pari_sp av = avma;
    3582             :   long nbgen, l1,l2;
    3583             :   GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
    3584         469 :   GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
    3585             : 
    3586         469 :   I1 = bid_get_ideal(bid1);
    3587         469 :   I2 = bid_get_ideal(bid2);
    3588         469 :   if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
    3589         259 :   G1 = bid_get_grp(bid1);
    3590         259 :   G2 = bid_get_grp(bid2);
    3591         259 :   x = idealmul(nf, I1,I2);
    3592         259 :   fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
    3593         259 :   fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
    3594         259 :   sprk1 = bid_get_sprk(bid1);
    3595         259 :   sprk2 = bid_get_sprk(bid2);
    3596         259 :   sprk = shallowconcat(sprk1, sprk2);
    3597             : 
    3598         259 :   cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
    3599         259 :   cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
    3600         259 :   gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
    3601         259 :   nbgen = l1+l2-2;
    3602         259 :   cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
    3603         259 :   if (nbgen)
    3604             :   {
    3605         259 :     GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
    3606           0 :     U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
    3607         259 :               : ZM_ZMV_mul(vecslice(U, 1, l1-1),   U1);
    3608           0 :     U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
    3609         259 :               : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
    3610         259 :     U = shallowconcat(U1, U2);
    3611             :   }
    3612             :   else
    3613           0 :     U = const_vec(lg(sprk), cgetg(1,t_MAT));
    3614             : 
    3615         259 :   if (gen)
    3616             :   {
    3617         259 :     GEN uv = zkchinese1init2(nf, I2, I1, x);
    3618         518 :     gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
    3619         259 :                         zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
    3620             :   }
    3621         259 :   sarch = bid_get_sarch(bid1); /* trivial */
    3622         259 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3623         259 :   x = mkvec2(x, bid_get_arch(bid1));
    3624         259 :   y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
    3625         259 :   return gerepilecopy(av,y);
    3626             : }
    3627             : 
    3628             : typedef struct _ideal_data {
    3629             :   GEN nf, emb, L, pr, prL, sgnU, archp;
    3630             : } ideal_data;
    3631             : 
    3632             : /* z <- ( z | f(v[i])_{i=1..#v} ) */
    3633             : static void
    3634      758382 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
    3635             : {
    3636      758382 :   long i, nz, lv = lg(v);
    3637             :   GEN z, Z;
    3638      758382 :   if (lv == 1) return;
    3639      222950 :   z = *pz; nz = lg(z)-1;
    3640      222950 :   *pz = Z = cgetg(lv + nz, typ(z));
    3641      371667 :   for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
    3642      223330 :   Z += nz;
    3643      491985 :   for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
    3644             : }
    3645             : static GEN
    3646         469 : join_idealinit(ideal_data *D, GEN x)
    3647         469 : { return join_bid(D->nf, x, D->prL); }
    3648             : static GEN
    3649      268442 : join_ideal(ideal_data *D, GEN x)
    3650      268442 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
    3651             : static GEN
    3652         448 : join_unit(ideal_data *D, GEN x)
    3653             : {
    3654         448 :   GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3655         448 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3656         448 :   return mkvec2(bid, v);
    3657             : }
    3658             : 
    3659             : GEN
    3660          49 : log_prk_units_init(GEN bnf)
    3661             : {
    3662          49 :   GEN U = bnf_has_fu(bnf);
    3663          49 :   if (U) U = matalgtobasis(bnf_get_nf(bnf), U);
    3664          21 :   else if (!(U = bnf_compactfu_mat(bnf))) (void)bnf_build_units(bnf);
    3665          49 :   return mkvec2(bnf_get_tuU(bnf), U);
    3666             : }
    3667             : /*  flag & nf_GEN : generators, otherwise no
    3668             :  *  flag &2 : units, otherwise no
    3669             :  *  flag &4 : ideals in HNF, otherwise bid
    3670             :  *  flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
    3671             : static GEN
    3672       11333 : Ideallist(GEN bnf, ulong bound, long flag)
    3673             : {
    3674       11333 :   const long do_units = flag & 2, big_id = !(flag & 4), cond = flag & 8;
    3675       11333 :   const long istar_flag = (flag & nf_GEN) | nf_INIT;
    3676             :   pari_sp av;
    3677             :   long i, j;
    3678       11333 :   GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
    3679             :   forprime_t S;
    3680             :   ideal_data ID;
    3681             :   GEN (*join_z)(ideal_data*, GEN);
    3682             : 
    3683       11333 :   if (do_units)
    3684             :   {
    3685          21 :     bnf = checkbnf(bnf);
    3686          21 :     nf = bnf_get_nf(bnf);
    3687          21 :     join_z = &join_unit;
    3688             :   }
    3689             :   else
    3690             :   {
    3691       11312 :     nf = checknf(bnf);
    3692       11312 :     join_z = big_id? &join_idealinit: &join_ideal;
    3693             :   }
    3694       11333 :   if ((long)bound <= 0) return empty;
    3695       11333 :   id = matid(nf_get_degree(nf));
    3696       11333 :   if (big_id) id = Idealstar(nf,id, istar_flag);
    3697             : 
    3698             :   /* z[i] will contain all "objects" of norm i. Depending on flag, this means
    3699             :    * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
    3700             :    * in vectors, computed one primary component at a time; join_z
    3701             :    * reconstructs the global object */
    3702       11333 :   BOUND = utoipos(bound);
    3703       11333 :   z = const_vec(bound, empty);
    3704       11333 :   U = do_units? log_prk_units_init(bnf): NULL;
    3705       11333 :   gel(z,1) = mkvec(U? mkvec2(id, empty): id);
    3706       11333 :   ID.nf = nf;
    3707             : 
    3708       11333 :   p = cgetipos(3);
    3709       11333 :   u_forprime_init(&S, 2, bound);
    3710       11333 :   av = avma;
    3711       92869 :   while ((p[2] = u_forprime_next(&S)))
    3712             :   {
    3713       81610 :     if (DEBUGLEVEL>1) err_printf("%ld ",p[2]);
    3714       81610 :     fa = idealprimedec_limit_norm(nf, p, BOUND);
    3715      163046 :     for (j = 1; j < lg(fa); j++)
    3716             :     {
    3717       81510 :       GEN pr = gel(fa,j), z2 = leafcopy(z);
    3718       81512 :       ulong Q, q = upr_norm(pr);
    3719             :       long l;
    3720       81511 :       ID.pr = ID.prL = pr;
    3721       81511 :       if (cond && q == 2) { l = 2; Q = 4; } else { l = 1; Q = q; }
    3722      184468 :       for (; Q <= bound; l++, Q *= q) /* add pr^l */
    3723             :       {
    3724             :         ulong iQ;
    3725      103035 :         ID.L = utoipos(l);
    3726      103034 :         if (big_id) {
    3727         210 :           ID.prL = Idealstarprk(nf, pr, l, istar_flag);
    3728         210 :           if (U)
    3729         189 :             ID.emb = Q == 2? empty
    3730         189 :                            : log_prk_units(nf, U, gel(bid_get_sprk(ID.prL),1));
    3731             :         }
    3732      861343 :         for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
    3733      758386 :           concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
    3734             :       }
    3735             :     }
    3736       81536 :     if (gc_needed(av,1))
    3737             :     {
    3738          18 :       if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
    3739          18 :       z = gerepilecopy(av, z);
    3740             :     }
    3741             :   }
    3742       11333 :   return z;
    3743             : }
    3744             : GEN
    3745          63 : gideallist(GEN bnf, GEN B, long flag)
    3746             : {
    3747          63 :   pari_sp av = avma;
    3748          63 :   if (typ(B) != t_INT)
    3749             :   {
    3750           0 :     B = gfloor(B);
    3751           0 :     if (typ(B) != t_INT) pari_err_TYPE("ideallist", B);
    3752           0 :     if (signe(B) < 0) B = gen_0;
    3753             :   }
    3754          63 :   if (signe(B) < 0)
    3755             :   {
    3756          28 :     if (flag != 4) pari_err_IMPL("ideallist with bid for single norm");
    3757          28 :     return gerepilecopy(av, ideals_by_norm(checknf(bnf), absi(B)));
    3758             :   }
    3759          35 :   if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
    3760          35 :   return gerepilecopy(av, Ideallist(bnf, itou(B), flag));
    3761             : }
    3762             : GEN
    3763       11298 : ideallist0(GEN bnf, long bound, long flag)
    3764             : {
    3765       11298 :   pari_sp av = avma;
    3766       11298 :   if (flag < 0 || flag > 15) pari_err_FLAG("ideallist");
    3767       11298 :   return gerepilecopy(av, Ideallist(bnf, bound, flag));
    3768             : }
    3769             : GEN
    3770       10563 : ideallist(GEN bnf,long bound) { return ideallist0(bnf,bound,4); }
    3771             : 
    3772             : /* bid = for module m (without arch. part), arch = archimedean part.
    3773             :  * Output: bid for [m,arch] */
    3774             : static GEN
    3775          42 : join_bid_arch(GEN nf, GEN bid, GEN archp)
    3776             : {
    3777          42 :   pari_sp av = avma;
    3778             :   GEN G, U;
    3779          42 :   GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
    3780             : 
    3781          42 :   checkbid(bid);
    3782          42 :   G = bid_get_grp(bid);
    3783          42 :   x = bid_get_ideal(bid);
    3784          42 :   sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
    3785          42 :   sprk = bid_get_sprk(bid);
    3786             : 
    3787          42 :   gen = (lg(G)>3)? gel(G,3): NULL;
    3788          42 :   cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
    3789          42 :   cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
    3790          42 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3791          42 :   U = split_U(U, sprk);
    3792          42 :   y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
    3793          42 :   return gerepilecopy(av,y);
    3794             : }
    3795             : static GEN
    3796          42 : join_arch(ideal_data *D, GEN x) {
    3797          42 :   return join_bid_arch(D->nf, x, D->archp);
    3798             : }
    3799             : static GEN
    3800          14 : join_archunit(ideal_data *D, GEN x) {
    3801          14 :   GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3802          14 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3803          14 :   return mkvec2(bid, v);
    3804             : }
    3805             : 
    3806             : /* L from ideallist, add archimedean part */
    3807             : GEN
    3808          14 : ideallistarch(GEN bnf, GEN L, GEN arch)
    3809             : {
    3810             :   pari_sp av;
    3811          14 :   long i, j, l = lg(L), lz;
    3812             :   GEN v, z, V, nf;
    3813             :   ideal_data ID;
    3814             :   GEN (*join_z)(ideal_data*, GEN);
    3815             : 
    3816          14 :   if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
    3817          14 :   if (l == 1) return cgetg(1,t_VEC);
    3818          14 :   z = gel(L,1);
    3819          14 :   if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3820          14 :   z = gel(z,1); /* either a bid or [bid,U] */
    3821          14 :   ID.archp = vec01_to_indices(arch);
    3822          14 :   if (lg(z) == 3)
    3823             :   { /* [bid,U]: do units */
    3824           7 :     bnf = checkbnf(bnf); nf = bnf_get_nf(bnf);
    3825           7 :     if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3826           7 :     ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
    3827           7 :     join_z = &join_archunit;
    3828             :   }
    3829             :   else
    3830             :   {
    3831           7 :     join_z = &join_arch;
    3832           7 :     nf = checknf(bnf);
    3833             :   }
    3834          14 :   ID.nf = nf;
    3835          14 :   av = avma; V = cgetg(l, t_VEC);
    3836          63 :   for (i = 1; i < l; i++)
    3837             :   {
    3838          49 :     z = gel(L,i); lz = lg(z);
    3839          49 :     gel(V,i) = v = cgetg(lz,t_VEC);
    3840          91 :     for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
    3841             :   }
    3842          14 :   return gerepilecopy(av,V);
    3843             : }

Generated by: LCOV version 1.16