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.12.1 lcov report (development 24988-2584e74448) Lines: 1727 1846 93.6 %
Date: 2020-01-26 05:57:03 Functions: 196 206 95.1 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8             : ANY WARRANTY WHATSOEVER.
       9             : 
      10             : Check the License for details. You should have received a copy of it, along
      11             : with the package; see the file 'COPYING'. If not, write to the Free Software
      12             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13             : 
      14             : /*******************************************************************/
      15             : /*                                                                 */
      16             : /*                       BASIC NF OPERATIONS                       */
      17             : /*                                                                 */
      18             : /*******************************************************************/
      19             : #include "pari.h"
      20             : #include "paripriv.h"
      21             : 
      22             : /*******************************************************************/
      23             : /*                                                                 */
      24             : /*                OPERATIONS OVER NUMBER FIELD ELEMENTS.           */
      25             : /*     represented as column vectors over the integral basis       */
      26             : /*                                                                 */
      27             : /*******************************************************************/
      28             : static GEN
      29    14363963 : get_tab(GEN nf, long *N)
      30             : {
      31    14363963 :   GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
      32    14363963 :   *N = nbrows(tab); return tab;
      33             : }
      34             : 
      35             : /* x != 0, y t_INT. Return x * y (not memory clean if x = 1) */
      36             : static GEN
      37   608544683 : _mulii(GEN x, GEN y) {
      38  1007766274 :   return is_pm1(x)? (signe(x) < 0)? negi(y): y
      39  1007770572 :                   : mulii(x, y);
      40             : }
      41             : 
      42             : GEN
      43       17864 : tablemul_ei_ej(GEN M, long i, long j)
      44             : {
      45             :   long N;
      46       17864 :   GEN tab = get_tab(M, &N);
      47       17864 :   tab += (i-1)*N; return gel(tab,j);
      48             : }
      49             : 
      50             : /* Outputs x.ei, where ei is the i-th elt of the algebra basis.
      51             :  * x an RgV of correct length and arbitrary content (polynomials, scalars...).
      52             :  * M is the multiplication table ei ej = sum_k M_k^(i,j) ek */
      53             : GEN
      54       10731 : tablemul_ei(GEN M, GEN x, long i)
      55             : {
      56             :   long j, k, N;
      57             :   GEN v, tab;
      58             : 
      59       10731 :   if (i==1) return gcopy(x);
      60       10731 :   tab = get_tab(M, &N);
      61       10731 :   if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
      62       10731 :   tab += (i-1)*N; v = cgetg(N+1,t_COL);
      63             :   /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
      64       74725 :   for (k=1; k<=N; k++)
      65             :   {
      66       63994 :     pari_sp av = avma;
      67       63994 :     GEN s = gen_0;
      68      458388 :     for (j=1; j<=N; j++)
      69             :     {
      70      394394 :       GEN c = gcoeff(tab,k,j);
      71      394394 :       if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
      72             :     }
      73       63994 :     gel(v,k) = gerepileupto(av,s);
      74             :   }
      75       10731 :   return v;
      76             : }
      77             : /* as tablemul_ei, assume x a ZV of correct length */
      78             : GEN
      79    10512018 : zk_ei_mul(GEN nf, GEN x, long i)
      80             : {
      81             :   long j, k, N;
      82             :   GEN v, tab;
      83             : 
      84    10512018 :   if (i==1) return ZC_copy(x);
      85    10512004 :   tab = get_tab(nf, &N); tab += (i-1)*N;
      86    10511987 :   v = cgetg(N+1,t_COL);
      87    73873182 :   for (k=1; k<=N; k++)
      88             :   {
      89    63361256 :     pari_sp av = avma;
      90    63361256 :     GEN s = gen_0;
      91   784908168 :     for (j=1; j<=N; j++)
      92             :     {
      93   721550097 :       GEN c = gcoeff(tab,k,j);
      94   721550097 :       if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
      95             :     }
      96    63358071 :     gel(v,k) = gerepileuptoint(av, s);
      97             :   }
      98    10511926 :   return v;
      99             : }
     100             : 
     101             : /* table of multiplication by wi in R[w1,..., wN] */
     102             : GEN
     103        2646 : ei_multable(GEN TAB, long i)
     104             : {
     105             :   long k,N;
     106        2646 :   GEN m, tab = get_tab(TAB, &N);
     107        2646 :   tab += (i-1)*N;
     108        2646 :   m = cgetg(N+1,t_MAT);
     109        2646 :   for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
     110        2646 :   return m;
     111             : }
     112             : 
     113             : GEN
     114     5297114 : zk_multable(GEN nf, GEN x)
     115             : {
     116     5297114 :   long i, l = lg(x);
     117     5297114 :   GEN mul = cgetg(l,t_MAT);
     118     5297142 :   gel(mul,1) = x; /* assume w_1 = 1 */
     119     5297142 :   for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
     120     5297060 :   return mul;
     121             : }
     122             : GEN
     123        1995 : multable(GEN M, GEN x)
     124             : {
     125             :   long i, N;
     126             :   GEN mul;
     127        1995 :   if (typ(x) == t_MAT) return x;
     128           0 :   M = get_tab(M, &N);
     129           0 :   if (typ(x) != t_COL) return scalarmat(x, N);
     130           0 :   mul = cgetg(N+1,t_MAT);
     131           0 :   gel(mul,1) = x; /* assume w_1 = 1 */
     132           0 :   for (i=2; i<=N; i++) gel(mul,i) = tablemul_ei(M,x,i);
     133           0 :   return mul;
     134             : }
     135             : 
     136             : /* x integral in nf; table of multiplication by x in ZK = Z[w1,..., wN].
     137             :  * Return a t_INT if x is scalar, and a ZM otherwise */
     138             : GEN
     139     3738274 : zk_scalar_or_multable(GEN nf, GEN x)
     140             : {
     141     3738274 :   long tx = typ(x);
     142     3738274 :   if (tx == t_MAT || tx == t_INT) return x;
     143     3711040 :   x = nf_to_scalar_or_basis(nf, x);
     144     3711040 :   return (typ(x) == t_COL)? zk_multable(nf, x): x;
     145             : }
     146             : 
     147             : GEN
     148       15330 : nftrace(GEN nf, GEN x)
     149             : {
     150       15330 :   pari_sp av = avma;
     151       15330 :   nf = checknf(nf);
     152       15330 :   x = nf_to_scalar_or_basis(nf, x);
     153       45969 :   x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
     154       30639 :                        : gmulgs(x, nf_get_degree(nf));
     155       15330 :   return gerepileupto(av, x);
     156             : }
     157             : GEN
     158         784 : rnfelttrace(GEN rnf, GEN x)
     159             : {
     160         784 :   pari_sp av = avma;
     161         784 :   checkrnf(rnf);
     162         784 :   x = rnfeltabstorel(rnf, x);
     163        2002 :   x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
     164        1309 :                           : gmulgs(x, rnf_get_degree(rnf));
     165         693 :   return gerepileupto(av, x);
     166             : }
     167             : 
     168             : /* assume nf is a genuine nf, fa a famat */
     169             : static GEN
     170           7 : famat_norm(GEN nf, GEN fa)
     171             : {
     172           7 :   pari_sp av = avma;
     173           7 :   GEN g = gel(fa,1), e = gel(fa,2), N = gen_1;
     174           7 :   long i, l = lg(g);
     175          21 :   for (i = 1; i < l; i++)
     176          14 :     N = gmul(N, powgi(nfnorm(nf, gel(g,i)), gel(e,i)));
     177           7 :   return gerepileupto(av, N);
     178             : }
     179             : GEN
     180       39889 : nfnorm(GEN nf, GEN x)
     181             : {
     182       39889 :   pari_sp av = avma;
     183       39889 :   nf = checknf(nf);
     184       39889 :   if (typ(x) == t_MAT) return famat_norm(nf, x);
     185       39882 :   x = nf_to_scalar_or_alg(nf, x);
     186      108782 :   x = (typ(x) == t_POL)? RgXQ_norm(x, nf_get_pol(nf))
     187       68900 :                        : gpowgs(x, nf_get_degree(nf));
     188       39882 :   return gerepileupto(av, x);
     189             : }
     190             : 
     191             : GEN
     192         231 : rnfeltnorm(GEN rnf, GEN x)
     193             : {
     194         231 :   pari_sp av = avma;
     195         231 :   checkrnf(rnf);
     196         231 :   x = rnfeltabstorel(rnf, x);
     197         378 :   x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gnorm(x))
     198         238 :                           : gpowgs(x, rnf_get_degree(rnf));
     199         140 :   return gerepileupto(av, x);
     200             : }
     201             : 
     202             : /* x + y in nf */
     203             : GEN
     204    16711198 : nfadd(GEN nf, GEN x, GEN y)
     205             : {
     206    16711198 :   pari_sp av = avma;
     207             :   GEN z;
     208             : 
     209    16711198 :   nf = checknf(nf);
     210    16711198 :   x = nf_to_scalar_or_basis(nf, x);
     211    16711198 :   y = nf_to_scalar_or_basis(nf, y);
     212    16711198 :   if (typ(x) != t_COL)
     213    13388515 :   { z = (typ(y) == t_COL)? RgC_Rg_add(y, x): gadd(x,y); }
     214             :   else
     215     3322683 :   { z = (typ(y) == t_COL)? RgC_add(x, y): RgC_Rg_add(x, y); }
     216    16711198 :   return gerepileupto(av, z);
     217             : }
     218             : /* x - y in nf */
     219             : GEN
     220     1249073 : nfsub(GEN nf, GEN x, GEN y)
     221             : {
     222     1249073 :   pari_sp av = avma;
     223             :   GEN z;
     224             : 
     225     1249073 :   nf = checknf(nf);
     226     1249073 :   x = nf_to_scalar_or_basis(nf, x);
     227     1249073 :   y = nf_to_scalar_or_basis(nf, y);
     228     1249073 :   if (typ(x) != t_COL)
     229      894950 :   { z = (typ(y) == t_COL)? Rg_RgC_sub(x,y): gsub(x,y); }
     230             :   else
     231      354123 :   { z = (typ(y) == t_COL)? RgC_sub(x,y): RgC_Rg_sub(x,y); }
     232     1249073 :   return gerepileupto(av, z);
     233             : }
     234             : 
     235             : /* product of ZC x,y in nf; ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
     236             : static GEN
     237     2215271 : nfmuli_ZC(GEN nf, GEN x, GEN y)
     238             : {
     239             :   long i, j, k, N;
     240     2215271 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     241             : 
     242    12166734 :   for (k = 1; k <= N; k++)
     243             :   {
     244     9951463 :     pari_sp av = avma;
     245     9951463 :     GEN s, TABi = TAB;
     246     9951463 :     if (k == 1)
     247     2215271 :       s = mulii(gel(x,1),gel(y,1));
     248             :     else
     249    15472384 :       s = addii(mulii(gel(x,1),gel(y,k)),
     250    15472384 :                 mulii(gel(x,k),gel(y,1)));
     251    83151975 :     for (i=2; i<=N; i++)
     252             :     {
     253    73200512 :       GEN t, xi = gel(x,i);
     254    73200512 :       TABi += N;
     255    73200512 :       if (!signe(xi)) continue;
     256             : 
     257    32254171 :       t = NULL;
     258   382860735 :       for (j=2; j<=N; j++)
     259             :       {
     260   350606564 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     261   350606564 :         if (!signe(c)) continue;
     262   126813106 :         p1 = _mulii(c, gel(y,j));
     263   126813106 :         t = t? addii(t, p1): p1;
     264             :       }
     265    32254171 :       if (t) s = addii(s, mulii(xi, t));
     266             :     }
     267     9951463 :     gel(v,k) = gerepileuptoint(av,s);
     268             :   }
     269     2215271 :   return v;
     270             : }
     271             : static int
     272    43141022 : is_famat(GEN x) { return typ(x) == t_MAT && lg(x) == 3; }
     273             : /* product of x and y in nf */
     274             : GEN
     275    22028532 : nfmul(GEN nf, GEN x, GEN y)
     276             : {
     277             :   GEN z;
     278    22028532 :   pari_sp av = avma;
     279             : 
     280    22028532 :   if (x == y) return nfsqr(nf,x);
     281             : 
     282    18849216 :   nf = checknf(nf);
     283    18849216 :   if (is_famat(x) || is_famat(y)) return famat_mul(x, y);
     284    18849125 :   x = nf_to_scalar_or_basis(nf, x);
     285    18849125 :   y = nf_to_scalar_or_basis(nf, y);
     286    18849125 :   if (typ(x) != t_COL)
     287             :   {
     288    14295269 :     if (isintzero(x)) return gen_0;
     289    10034950 :     z = (typ(y) == t_COL)? RgC_Rg_mul(y, x): gmul(x,y); }
     290             :   else
     291             :   {
     292     4553856 :     if (typ(y) != t_COL)
     293             :     {
     294     3186477 :       if (isintzero(y)) return gen_0;
     295      700896 :       z = RgC_Rg_mul(x, y);
     296             :     }
     297             :     else
     298             :     {
     299             :       GEN dx, dy;
     300     1367379 :       x = Q_remove_denom(x, &dx);
     301     1367379 :       y = Q_remove_denom(y, &dy);
     302     1367379 :       z = nfmuli_ZC(nf,x,y);
     303     1367379 :       dx = mul_denom(dx,dy);
     304     1367379 :       if (dx) z = ZC_Z_div(z, dx);
     305             :     }
     306             :   }
     307    12103225 :   return gerepileupto(av, z);
     308             : }
     309             : /* square of ZC x in nf */
     310             : static GEN
     311     1605457 : nfsqri_ZC(GEN nf, GEN x)
     312             : {
     313             :   long i, j, k, N;
     314     1605457 :   GEN TAB = get_tab(nf, &N), v = cgetg(N+1,t_COL);
     315             : 
     316    13801575 :   for (k = 1; k <= N; k++)
     317             :   {
     318    12196118 :     pari_sp av = avma;
     319    12196118 :     GEN s, TABi = TAB;
     320    12196118 :     if (k == 1)
     321     1605457 :       s = sqri(gel(x,1));
     322             :     else
     323    10590661 :       s = shifti(mulii(gel(x,1),gel(x,k)), 1);
     324   141884588 :     for (i=2; i<=N; i++)
     325             :     {
     326   129688470 :       GEN p1, c, t, xi = gel(x,i);
     327   129688470 :       TABi += N;
     328   129688470 :       if (!signe(xi)) continue;
     329             : 
     330    43354980 :       c = gcoeff(TABi, k, i);
     331    43354980 :       t = signe(c)? _mulii(c,xi): NULL;
     332   480555420 :       for (j=i+1; j<=N; j++)
     333             :       {
     334   437200440 :         c = gcoeff(TABi, k, j);
     335   437200440 :         if (!signe(c)) continue;
     336   184664011 :         p1 = _mulii(c, shifti(gel(x,j),1));
     337   184664011 :         t = t? addii(t, p1): p1;
     338             :       }
     339    43354980 :       if (t) s = addii(s, mulii(xi, t));
     340             :     }
     341    12196118 :     gel(v,k) = gerepileuptoint(av,s);
     342             :   }
     343     1605457 :   return v;
     344             : }
     345             : /* square of x in nf */
     346             : GEN
     347     5053805 : nfsqr(GEN nf, GEN x)
     348             : {
     349     5053805 :   pari_sp av = avma;
     350             :   GEN z;
     351             : 
     352     5053805 :   nf = checknf(nf);
     353     5053805 :   if (is_famat(x)) return famat_sqr(x);
     354     5053805 :   x = nf_to_scalar_or_basis(nf, x);
     355     5053805 :   if (typ(x) != t_COL) z = gsqr(x);
     356             :   else
     357             :   {
     358             :     GEN dx;
     359       95964 :     x = Q_remove_denom(x, &dx);
     360       95964 :     z = nfsqri_ZC(nf,x);
     361       95964 :     if (dx) z = RgC_Rg_div(z, sqri(dx));
     362             :   }
     363     5053805 :   return gerepileupto(av, z);
     364             : }
     365             : 
     366             : /* x a ZC, v a t_COL of ZC/Z */
     367             : GEN
     368      136188 : zkC_multable_mul(GEN v, GEN x)
     369             : {
     370      136188 :   long i, l = lg(v);
     371      136188 :   GEN y = cgetg(l, t_COL);
     372      477340 :   for (i = 1; i < l; i++)
     373             :   {
     374      341152 :     GEN c = gel(v,i);
     375      341152 :     if (typ(c)!=t_COL) {
     376           0 :       if (!isintzero(c)) c = ZC_Z_mul(gel(x,1), c);
     377             :     } else {
     378      341152 :       c = ZM_ZC_mul(x,c);
     379      341152 :       if (ZV_isscalar(c)) c = gel(c,1);
     380             :     }
     381      341152 :     gel(y,i) = c;
     382             :   }
     383      136188 :   return y;
     384             : }
     385             : 
     386             : GEN
     387       56133 : nfC_multable_mul(GEN v, GEN x)
     388             : {
     389       56133 :   long i, l = lg(v);
     390       56133 :   GEN y = cgetg(l, t_COL);
     391      364000 :   for (i = 1; i < l; i++)
     392             :   {
     393      307867 :     GEN c = gel(v,i);
     394      307867 :     if (typ(c)!=t_COL) {
     395      249382 :       if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
     396             :     } else {
     397       58485 :       c = RgM_RgC_mul(x,c);
     398       58485 :       if (QV_isscalar(c)) c = gel(c,1);
     399             :     }
     400      307867 :     gel(y,i) = c;
     401             :   }
     402       56133 :   return y;
     403             : }
     404             : 
     405             : GEN
     406      183610 : nfC_nf_mul(GEN nf, GEN v, GEN x)
     407             : {
     408             :   long tx;
     409             :   GEN y;
     410             : 
     411      183610 :   x = nf_to_scalar_or_basis(nf, x);
     412      183610 :   tx = typ(x);
     413      183610 :   if (tx != t_COL)
     414             :   {
     415             :     long l, i;
     416      134421 :     if (tx == t_INT)
     417             :     {
     418      125629 :       long s = signe(x);
     419      125629 :       if (!s) return zerocol(lg(v)-1);
     420      118776 :       if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
     421             :     }
     422       43701 :     l = lg(v); y = cgetg(l, t_COL);
     423      313313 :     for (i=1; i < l; i++)
     424             :     {
     425      269612 :       GEN c = gel(v,i);
     426      269612 :       if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
     427      269612 :       gel(y,i) = c;
     428             :     }
     429       43701 :     return y;
     430             :   }
     431             :   else
     432             :   {
     433             :     GEN dx;
     434       49189 :     x = zk_multable(nf, Q_remove_denom(x,&dx));
     435       49189 :     y = nfC_multable_mul(v, x);
     436       49189 :     return dx? RgC_Rg_div(y, dx): y;
     437             :   }
     438             : }
     439             : static GEN
     440        9079 : mulbytab(GEN M, GEN c)
     441        9079 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
     442             : GEN
     443        1995 : tablemulvec(GEN M, GEN x, GEN v)
     444             : {
     445             :   long l, i;
     446             :   GEN y;
     447             : 
     448        1995 :   if (typ(x) == t_COL && RgV_isscalar(x))
     449             :   {
     450           0 :     x = gel(x,1);
     451           0 :     return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
     452             :   }
     453        1995 :   x = multable(M, x); /* multiplication table by x */
     454        1995 :   y = cgetg_copy(v, &l);
     455        1995 :   if (typ(v) == t_POL)
     456             :   {
     457        1995 :     y[1] = v[1];
     458        1995 :     for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     459        1995 :     y = normalizepol(y);
     460             :   }
     461             :   else
     462             :   {
     463           0 :     for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     464             :   }
     465        1995 :   return y;
     466             : }
     467             : 
     468             : GEN
     469      371484 : zkmultable_capZ(GEN mx) { return Q_denom(zkmultable_inv(mx)); }
     470             : GEN
     471      432634 : zkmultable_inv(GEN mx) { return ZM_gauss(mx, col_ei(lg(mx)-1,1)); }
     472             : /* nf a true nf, x a ZC */
     473             : GEN
     474       61150 : zk_inv(GEN nf, GEN x) { return zkmultable_inv(zk_multable(nf,x)); }
     475             : 
     476             : /* inverse of x in nf */
     477             : GEN
     478       69643 : nfinv(GEN nf, GEN x)
     479             : {
     480       69643 :   pari_sp av = avma;
     481             :   GEN z;
     482             : 
     483       69643 :   nf = checknf(nf);
     484       69643 :   if (is_famat(x)) return famat_inv(x);
     485       69643 :   x = nf_to_scalar_or_basis(nf, x);
     486       69643 :   if (typ(x) == t_COL)
     487             :   {
     488             :     GEN d;
     489       28091 :     x = Q_remove_denom(x, &d);
     490       28091 :     z = zk_inv(nf, x);
     491       28091 :     if (d) z = RgC_Rg_mul(z, d);
     492             :   }
     493             :   else
     494       41552 :     z = ginv(x);
     495       69643 :   return gerepileupto(av, z);
     496             : }
     497             : 
     498             : /* quotient of x and y in nf */
     499             : GEN
     500       33489 : nfdiv(GEN nf, GEN x, GEN y)
     501             : {
     502       33489 :   pari_sp av = avma;
     503             :   GEN z;
     504             : 
     505       33489 :   nf = checknf(nf);
     506       33489 :   if (is_famat(x) || is_famat(y)) return famat_div(x,y);
     507       33398 :   y = nf_to_scalar_or_basis(nf, y);
     508       33398 :   if (typ(y) != t_COL)
     509             :   {
     510       20258 :     x = nf_to_scalar_or_basis(nf, x);
     511       20258 :     z = (typ(x) == t_COL)? RgC_Rg_div(x, y): gdiv(x,y);
     512             :   }
     513             :   else
     514             :   {
     515             :     GEN d;
     516       13140 :     y = Q_remove_denom(y, &d);
     517       13140 :     z = nfmul(nf, x, zk_inv(nf,y));
     518       13140 :     if (d) z = RgC_Rg_mul(z, d);
     519             :   }
     520       33398 :   return gerepileupto(av, z);
     521             : }
     522             : 
     523             : /* product of INTEGERS (t_INT or ZC) x and y in nf */
     524             : GEN
     525     1076490 : nfmuli(GEN nf, GEN x, GEN y)
     526             : {
     527     1076490 :   if (typ(x) == t_INT) return (typ(y) == t_COL)? ZC_Z_mul(y, x): mulii(x,y);
     528      896223 :   if (typ(y) == t_INT) return ZC_Z_mul(x, y);
     529      847892 :   return nfmuli_ZC(nf, x, y);
     530             : }
     531             : GEN
     532     1509493 : nfsqri(GEN nf, GEN x)
     533     1509493 : { return (typ(x) == t_INT)? sqri(x): nfsqri_ZC(nf, x); }
     534             : 
     535             : /* both x and y are RgV */
     536             : GEN
     537           0 : tablemul(GEN TAB, GEN x, GEN y)
     538             : {
     539             :   long i, j, k, N;
     540             :   GEN s, v;
     541           0 :   if (typ(x) != t_COL) return gmul(x, y);
     542           0 :   if (typ(y) != t_COL) return gmul(y, x);
     543           0 :   N = lg(x)-1;
     544           0 :   v = cgetg(N+1,t_COL);
     545           0 :   for (k=1; k<=N; k++)
     546             :   {
     547           0 :     pari_sp av = avma;
     548           0 :     GEN TABi = TAB;
     549           0 :     if (k == 1)
     550           0 :       s = gmul(gel(x,1),gel(y,1));
     551             :     else
     552           0 :       s = gadd(gmul(gel(x,1),gel(y,k)),
     553           0 :                gmul(gel(x,k),gel(y,1)));
     554           0 :     for (i=2; i<=N; i++)
     555             :     {
     556           0 :       GEN t, xi = gel(x,i);
     557           0 :       TABi += N;
     558           0 :       if (gequal0(xi)) continue;
     559             : 
     560           0 :       t = NULL;
     561           0 :       for (j=2; j<=N; j++)
     562             :       {
     563           0 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     564           0 :         if (gequal0(c)) continue;
     565           0 :         p1 = gmul(c, gel(y,j));
     566           0 :         t = t? gadd(t, p1): p1;
     567             :       }
     568           0 :       if (t) s = gadd(s, gmul(xi, t));
     569             :     }
     570           0 :     gel(v,k) = gerepileupto(av,s);
     571             :   }
     572           0 :   return v;
     573             : }
     574             : GEN
     575       42091 : tablesqr(GEN TAB, GEN x)
     576             : {
     577             :   long i, j, k, N;
     578             :   GEN s, v;
     579             : 
     580       42091 :   if (typ(x) != t_COL) return gsqr(x);
     581       42091 :   N = lg(x)-1;
     582       42091 :   v = cgetg(N+1,t_COL);
     583             : 
     584      298795 :   for (k=1; k<=N; k++)
     585             :   {
     586      256704 :     pari_sp av = avma;
     587      256704 :     GEN TABi = TAB;
     588      256704 :     if (k == 1)
     589       42091 :       s = gsqr(gel(x,1));
     590             :     else
     591      214613 :       s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
     592     1611036 :     for (i=2; i<=N; i++)
     593             :     {
     594     1354332 :       GEN p1, c, t, xi = gel(x,i);
     595     1354332 :       TABi += N;
     596     1354332 :       if (gequal0(xi)) continue;
     597             : 
     598      372918 :       c = gcoeff(TABi, k, i);
     599      372918 :       t = !gequal0(c)? gmul(c,xi): NULL;
     600     1497363 :       for (j=i+1; j<=N; j++)
     601             :       {
     602     1124445 :         c = gcoeff(TABi, k, j);
     603     1124445 :         if (gequal0(c)) continue;
     604      586551 :         p1 = gmul(gmul2n(c,1), gel(x,j));
     605      586551 :         t = t? gadd(t, p1): p1;
     606             :       }
     607      372918 :       if (t) s = gadd(s, gmul(xi, t));
     608             :     }
     609      256704 :     gel(v,k) = gerepileupto(av,s);
     610             :   }
     611       42091 :   return v;
     612             : }
     613             : 
     614             : static GEN
     615      163337 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
     616             : static GEN
     617      387242 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
     618             : 
     619             : /* Compute z^n in nf, left-shift binary powering */
     620             : GEN
     621      252108 : nfpow(GEN nf, GEN z, GEN n)
     622             : {
     623      252108 :   pari_sp av = avma;
     624             :   long s;
     625             :   GEN x, cx;
     626             : 
     627      252108 :   if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
     628      252108 :   nf = checknf(nf);
     629      252108 :   s = signe(n); if (!s) return gen_1;
     630      252108 :   if (is_famat(z)) return famat_pow(z, n);
     631      252094 :   x = nf_to_scalar_or_basis(nf, z);
     632      252094 :   if (typ(x) != t_COL) return powgi(x,n);
     633      247997 :   if (s < 0)
     634             :   { /* simplified nfinv */
     635             :     GEN d;
     636       13866 :     x = Q_remove_denom(x, &d);
     637       13866 :     x = zk_inv(nf, x);
     638       13866 :     x = primitive_part(x, &cx);
     639       13866 :     cx = mul_content(cx, d);
     640       13866 :     n = negi(n);
     641             :   }
     642             :   else
     643      234131 :     x = primitive_part(x, &cx);
     644      247997 :   x = gen_pow_i(x, n, (void*)nf, _sqr, _mul);
     645      247997 :   if (cx)
     646       53029 :     x = gerepileupto(av, gmul(x, powgi(cx, n)));
     647             :   else
     648      194968 :     x = gerepilecopy(av, x);
     649      247997 :   return x;
     650             : }
     651             : /* Compute z^n in nf, left-shift binary powering */
     652             : GEN
     653       55832 : nfpow_u(GEN nf, GEN z, ulong n)
     654             : {
     655       55832 :   pari_sp av = avma;
     656             :   GEN x, cx;
     657             : 
     658       55832 :   nf = checknf(nf);
     659       55832 :   if (!n) return gen_1;
     660       55832 :   x = nf_to_scalar_or_basis(nf, z);
     661       55832 :   if (typ(x) != t_COL) return gpowgs(x,n);
     662       27258 :   x = primitive_part(x, &cx);
     663       27258 :   x = gen_powu_i(x, n, (void*)nf, _sqr, _mul);
     664       27258 :   if (cx)
     665             :   {
     666        1113 :     x = gmul(x, powgi(cx, utoipos(n)));
     667        1113 :     return gerepileupto(av,x);
     668             :   }
     669       26145 :   return gerepilecopy(av, x);
     670             : }
     671             : 
     672             : static GEN
     673     2570687 : _nf_red(void *E, GEN x) { (void)E; return x; }
     674             : 
     675             : static GEN
     676    10367140 : _nf_add(void *E, GEN x, GEN y) { return nfadd((GEN)E,x,y); }
     677             : 
     678             : static GEN
     679      623763 : _nf_neg(void *E, GEN x) { (void)E; return gneg(x); }
     680             : 
     681             : static GEN
     682    12473391 : _nf_mul(void *E, GEN x, GEN y) { return nfmul((GEN)E,x,y); }
     683             : 
     684             : static GEN
     685       43218 : _nf_inv(void *E, GEN x) { return nfinv((GEN)E,x); }
     686             : 
     687             : static GEN
     688        8498 : _nf_s(void *E, long x) { (void)E; return stoi(x); }
     689             : 
     690             : static const struct bb_field nf_field={_nf_red,_nf_add,_nf_mul,_nf_neg,
     691             :                                         _nf_inv,&gequal0,_nf_s };
     692             : 
     693      186984 : const struct bb_field *get_nf_field(void **E, GEN nf)
     694      186984 : { *E = (void*)nf; return &nf_field; }
     695             : 
     696             : GEN
     697          14 : nfM_det(GEN nf, GEN M)
     698             : {
     699             :   void *E;
     700          14 :   const struct bb_field *S = get_nf_field(&E, nf);
     701          14 :   return gen_det(M, E, S);
     702             : }
     703             : GEN
     704        8484 : nfM_inv(GEN nf, GEN M)
     705             : {
     706             :   void *E;
     707        8484 :   const struct bb_field *S = get_nf_field(&E, nf);
     708        8484 :   return gen_Gauss(M, matid(lg(M)-1), E, S);
     709             : }
     710             : GEN
     711        8232 : nfM_mul(GEN nf, GEN A, GEN B)
     712             : {
     713             :   void *E;
     714        8232 :   const struct bb_field *S = get_nf_field(&E, nf);
     715        8232 :   return gen_matmul(A, B, E, S);
     716             : }
     717             : GEN
     718      170254 : nfM_nfC_mul(GEN nf, GEN A, GEN B)
     719             : {
     720             :   void *E;
     721      170254 :   const struct bb_field *S = get_nf_field(&E, nf);
     722      170254 :   return gen_matcolmul(A, B, E, S);
     723             : }
     724             : 
     725             : /* valuation of integral x (ZV), with resp. to prime ideal pr */
     726             : long
     727     9099369 : ZC_nfvalrem(GEN x, GEN pr, GEN *newx)
     728             : {
     729             :   long i, v, l;
     730     9099369 :   GEN r, y, p = pr_get_p(pr), mul = pr_get_tau(pr);
     731             : 
     732             :   /* p inert */
     733     9099374 :   if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
     734     9084884 :   y = cgetg_copy(x, &l); /* will hold the new x */
     735     9084895 :   x = leafcopy(x);
     736    13116200 :   for(v=0;; v++)
     737             :   {
     738    48965552 :     for (i=1; i<l; i++)
     739             :     { /* is (x.b)[i] divisible by p ? */
     740    40902956 :       gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
     741    40902930 :       if (r != gen_0) { if (newx) *newx = x; return v; }
     742             :     }
     743     4031298 :     swap(x, y);
     744             :   }
     745             : }
     746             : long
     747     8690581 : ZC_nfval(GEN x, GEN P)
     748     8690581 : { return ZC_nfvalrem(x, P, NULL); }
     749             : 
     750             : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
     751             : int
     752      385312 : ZC_prdvd(GEN x, GEN P)
     753             : {
     754      385312 :   pari_sp av = avma;
     755             :   long i, l;
     756      385312 :   GEN p = pr_get_p(P), mul = pr_get_tau(P);
     757      385312 :   if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
     758      384969 :   l = lg(x);
     759     1550795 :   for (i=1; i<l; i++)
     760     1375228 :     if (!dvdii(ZMrow_ZC_mul(mul,x,i), p)) return gc_bool(av,0);
     761      175567 :   return gc_bool(av,1);
     762             : }
     763             : 
     764             : int
     765          28 : pr_equal(GEN P, GEN Q)
     766             : {
     767          28 :   GEN gQ, p = pr_get_p(P);
     768          28 :   long e = pr_get_e(P), f = pr_get_f(P), n;
     769          28 :   if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
     770          14 :     return 0;
     771          14 :   gQ = pr_get_gen(Q); n = lg(gQ)-1;
     772          14 :   if (2*e*f > n) return 1; /* room for only one such pr */
     773           7 :   return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(gQ, P);
     774             : }
     775             : 
     776             : GEN
     777          98 : famat_nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     778             : {
     779          98 :   pari_sp av = avma;
     780          98 :   GEN P = gel(x,1), E = gel(x,2), V = gen_0, y = NULL;
     781          98 :   long l = lg(P), i;
     782          98 :   if (py) { *py = gen_1; y = cgetg(l, t_COL); }
     783             : 
     784         294 :   for (i = 1; i < l; i++)
     785             :   {
     786         196 :     GEN e = gel(E,i);
     787             :     long v;
     788         196 :     if (!signe(e)) continue;
     789         196 :     v = nfvalrem(nf, gel(P,i), pr, py? &gel(y,i): NULL);
     790         196 :     if (v == LONG_MAX) { set_avma(av); if (py) *py = gen_0; return mkoo(); }
     791         196 :     V = addmulii(V, stoi(v), e);
     792             :   }
     793          98 :   if (!py) V = gerepileuptoint(av, V);
     794             :   else
     795             :   {
     796          28 :     y = mkmat2(y, gel(x,2));
     797          28 :     gerepileall(av, 2, &V, &y); *py = y;
     798             :   }
     799          98 :   return V;
     800             : }
     801             : long
     802     1306591 : nfval(GEN nf, GEN x, GEN pr)
     803             : {
     804     1306591 :   pari_sp av = avma;
     805             :   long w, e;
     806             :   GEN cx, p;
     807             : 
     808     1306591 :   if (gequal0(x)) return LONG_MAX;
     809     1305004 :   nf = checknf(nf);
     810     1305004 :   checkprid(pr);
     811     1305001 :   p = pr_get_p(pr);
     812     1304997 :   e = pr_get_e(pr);
     813     1304998 :   x = nf_to_scalar_or_basis(nf, x);
     814     1304978 :   if (typ(x) != t_COL) return e*Q_pval(x,p);
     815      108745 :   x = Q_primitive_part(x, &cx);
     816      108731 :   w = ZC_nfval(x,pr);
     817      108736 :   if (cx) w += e*Q_pval(cx,p);
     818      108758 :   return gc_long(av,w);
     819             : }
     820             : 
     821             : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
     822             : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
     823             : static GEN
     824       20174 : powp(GEN nf, GEN pr, long v)
     825             : {
     826             :   GEN b, z;
     827             :   long e;
     828       20174 :   if (!v) return gen_1;
     829       20048 :   b = pr_get_tau(pr);
     830       20048 :   if (typ(b) == t_INT) return gen_1;
     831        1344 :   e = pr_get_e(pr);
     832        1344 :   z = gel(b,1);
     833        1344 :   if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
     834        1344 :   if (v < 0) { v = -v; z = nfinv(nf, z); }
     835        1344 :   if (v != 1) z = nfpow_u(nf, z, v);
     836        1344 :   return z;
     837             : }
     838             : long
     839       64869 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     840             : {
     841       64869 :   pari_sp av = avma;
     842             :   long w, e;
     843             :   GEN cx, p, t;
     844             : 
     845       64869 :   if (!py) return nfval(nf,x,pr);
     846       64610 :   if (gequal0(x)) { *py = gen_0; return LONG_MAX; }
     847       64554 :   nf = checknf(nf);
     848       64554 :   checkprid(pr);
     849       64554 :   p = pr_get_p(pr);
     850       64554 :   e = pr_get_e(pr);
     851       64554 :   x = nf_to_scalar_or_basis(nf, x);
     852       64554 :   if (typ(x) != t_COL) {
     853       52885 :     w = Q_pvalrem(x,p, py);
     854       52885 :     if (!w) { *py = gerepilecopy(av, x); return 0; }
     855       18984 :     *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
     856       18984 :     return e*w;
     857             :   }
     858       11669 :   x = Q_primitive_part(x, &cx);
     859       11669 :   w = ZC_nfvalrem(x,pr, py);
     860       11669 :   if (cx)
     861             :   {
     862        1190 :     long v = Q_pvalrem(cx,p, &t);
     863        1190 :     *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
     864        1190 :     *py = gerepileupto(av, *py);
     865        1190 :     w += e*v;
     866             :   }
     867             :   else
     868       10479 :     *py = gerepilecopy(av, *py);
     869       11669 :   return w;
     870             : }
     871             : GEN
     872         154 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     873             : {
     874             :   long v;
     875         154 :   if (is_famat(x)) return famat_nfvalrem(nf, x, pr, py);
     876         147 :   v = nfvalrem(nf,x,pr,py);
     877         147 :   return v == LONG_MAX? mkoo(): stoi(v);
     878             : }
     879             : 
     880             : /* true nf */
     881             : GEN
     882       98077 : coltoalg(GEN nf, GEN x)
     883             : {
     884       98077 :   return mkpolmod( nf_to_scalar_or_alg(nf, x), nf_get_pol(nf) );
     885             : }
     886             : 
     887             : GEN
     888      144361 : basistoalg(GEN nf, GEN x)
     889             : {
     890             :   GEN T;
     891             : 
     892      144361 :   nf = checknf(nf);
     893      144361 :   switch(typ(x))
     894             :   {
     895             :     case t_COL: {
     896       91938 :       pari_sp av = avma;
     897       91938 :       return gerepilecopy(av, coltoalg(nf, x));
     898             :     }
     899             :     case t_POLMOD:
     900       32039 :       T = nf_get_pol(nf);
     901       32039 :       if (!RgX_equal_var(T,gel(x,1)))
     902           0 :         pari_err_MODULUS("basistoalg", T,gel(x,1));
     903       32039 :       return gcopy(x);
     904             :     case t_POL:
     905        1848 :       T = nf_get_pol(nf);
     906        1848 :       if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
     907        1848 :       retmkpolmod(RgX_rem(x, T), ZX_copy(T));
     908             :     case t_INT:
     909             :     case t_FRAC:
     910       18536 :       T = nf_get_pol(nf);
     911       18536 :       retmkpolmod(gcopy(x), ZX_copy(T));
     912             :     default:
     913           0 :       pari_err_TYPE("basistoalg",x);
     914             :       return NULL; /* LCOV_EXCL_LINE */
     915             :   }
     916             : }
     917             : 
     918             : /* true nf, x a t_POL */
     919             : static GEN
     920     1596724 : pol_to_scalar_or_basis(GEN nf, GEN x)
     921             : {
     922     1596724 :   GEN T = nf_get_pol(nf);
     923     1596726 :   long l = lg(x);
     924     1596726 :   if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
     925     1596667 :   if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
     926     1596667 :   if (l == 2) return gen_0;
     927      961354 :   if (l == 3)
     928             :   {
     929      213518 :     x = gel(x,2);
     930      213518 :     if (!is_rational_t(typ(x))) pari_err_TYPE("nf_to_scalar_or_basis",x);
     931      213510 :     return x;
     932             :   }
     933      747836 :   return poltobasis(nf,x);
     934             : }
     935             : /* Assume nf is a genuine nf. */
     936             : GEN
     937    89942174 : nf_to_scalar_or_basis(GEN nf, GEN x)
     938             : {
     939    89942174 :   switch(typ(x))
     940             :   {
     941             :     case t_INT: case t_FRAC:
     942    66347206 :       return x;
     943             :     case t_POLMOD:
     944      193742 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
     945      193688 :       switch(typ(x))
     946             :       {
     947       41306 :         case t_INT: case t_FRAC: return x;
     948      152382 :         case t_POL: return pol_to_scalar_or_basis(nf,x);
     949             :       }
     950           0 :       break;
     951     1444358 :     case t_POL: return pol_to_scalar_or_basis(nf,x);
     952             :     case t_COL:
     953    21956899 :       if (lg(x)-1 != nf_get_degree(nf)) break;
     954    21956835 :       return QV_isscalar(x)? gel(x,1): x;
     955             :   }
     956          54 :   pari_err_TYPE("nf_to_scalar_or_basis",x);
     957             :   return NULL; /* LCOV_EXCL_LINE */
     958             : }
     959             : /* Let x be a polynomial with coefficients in Q or nf. Return the same
     960             :  * polynomial with coefficients expressed as vectors (on the integral basis).
     961             :  * No consistency checks, not memory-clean. */
     962             : GEN
     963        7154 : RgX_to_nfX(GEN nf, GEN x)
     964             : {
     965             :   long i, l;
     966        7154 :   GEN y = cgetg_copy(x, &l); y[1] = x[1];
     967        7154 :   for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
     968        7154 :   return y;
     969             : }
     970             : 
     971             : /* Assume nf is a genuine nf. */
     972             : GEN
     973     1163342 : nf_to_scalar_or_alg(GEN nf, GEN x)
     974             : {
     975     1163342 :   switch(typ(x))
     976             :   {
     977             :     case t_INT: case t_FRAC:
     978       37371 :       return x;
     979             :     case t_POLMOD:
     980         455 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
     981         455 :       if (typ(x) != t_POL) return x;
     982             :       /* fall through */
     983             :     case t_POL:
     984             :     {
     985       20812 :       GEN T = nf_get_pol(nf);
     986       20812 :       long l = lg(x);
     987       20812 :       if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
     988       20812 :       if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
     989       20812 :       if (l == 2) return gen_0;
     990       20812 :       if (l == 3) return gel(x,2);
     991       20784 :       return x;
     992             :     }
     993             :     case t_COL:
     994             :     {
     995             :       GEN dx;
     996     1105104 :       if (lg(x)-1 != nf_get_degree(nf)) break;
     997     2210222 :       if (QV_isscalar(x)) return gel(x,1);
     998     1068709 :       x = Q_remove_denom(x, &dx);
     999     1068709 :       x = RgV_RgC_mul(nf_get_zkprimpart(nf), x);
    1000     1068720 :       dx = mul_denom(dx, nf_get_zkden(nf));
    1001     1068714 :       return gdiv(x,dx);
    1002             :     }
    1003             :   }
    1004          48 :   pari_err_TYPE("nf_to_scalar_or_alg",x);
    1005             :   return NULL; /* LCOV_EXCL_LINE */
    1006             : }
    1007             : 
    1008             : /* gmul(A, RgX_to_RgC(x)), A t_MAT of compatible dimensions */
    1009             : GEN
    1010        1337 : RgM_RgX_mul(GEN A, GEN x)
    1011             : {
    1012        1337 :   long i, l = lg(x)-1;
    1013             :   GEN z;
    1014        1337 :   if (l == 1) return zerocol(nbrows(A));
    1015        1337 :   z = gmul(gel(x,2), gel(A,1));
    1016        2541 :   for (i = 2; i < l; i++)
    1017        1204 :     if (!gequal0(gel(x,i+1))) z = gadd(z, gmul(gel(x,i+1), gel(A,i)));
    1018        1337 :   return z;
    1019             : }
    1020             : GEN
    1021     2981993 : ZM_ZX_mul(GEN A, GEN x)
    1022             : {
    1023     2981993 :   long i, l = lg(x)-1;
    1024             :   GEN z;
    1025     2981993 :   if (l == 1) return zerocol(nbrows(A));
    1026     2981107 :   z = ZC_Z_mul(gel(A,1), gel(x,2));
    1027    11911031 :   for (i = 2; i < l ; i++)
    1028     8930251 :     if (signe(gel(x,i+1))) z = ZC_add(z, ZC_Z_mul(gel(A,i), gel(x,i+1)));
    1029     2980780 :   return z;
    1030             : }
    1031             : /* x a t_POL, nf a genuine nf. No garbage collecting. No check.  */
    1032             : GEN
    1033     2758204 : poltobasis(GEN nf, GEN x)
    1034             : {
    1035     2758204 :   GEN d, T = nf_get_pol(nf);
    1036     2758174 :   if (varn(x) != varn(T)) pari_err_VAR( "poltobasis", x,T);
    1037     2758041 :   if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1038     2757994 :   x = Q_remove_denom(x, &d);
    1039     2757953 :   if (!RgX_is_ZX(x)) pari_err_TYPE("poltobasis",x);
    1040     2757928 :   x = ZM_ZX_mul(nf_get_invzk(nf), x);
    1041     2757579 :   if (d) x = RgC_Rg_div(x, d);
    1042     2757587 :   return x;
    1043             : }
    1044             : 
    1045             : GEN
    1046      300601 : algtobasis(GEN nf, GEN x)
    1047             : {
    1048             :   pari_sp av;
    1049             : 
    1050      300601 :   nf = checknf(nf);
    1051      300601 :   switch(typ(x))
    1052             :   {
    1053             :     case t_POLMOD:
    1054      107107 :       if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
    1055           7 :         pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
    1056      107100 :       x = gel(x,2);
    1057      107100 :       switch(typ(x))
    1058             :       {
    1059             :         case t_INT:
    1060        7497 :         case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1061             :         case t_POL:
    1062       99603 :           av = avma;
    1063       99603 :           return gerepileupto(av,poltobasis(nf,x));
    1064             :       }
    1065           0 :       break;
    1066             : 
    1067             :     case t_POL:
    1068       99337 :       av = avma;
    1069       99337 :       return gerepileupto(av,poltobasis(nf,x));
    1070             : 
    1071             :     case t_COL:
    1072       21693 :       if (!RgV_is_QV(x)) pari_err_TYPE("nfalgtobasis",x);
    1073       21686 :       if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
    1074       21686 :       return gcopy(x);
    1075             : 
    1076             :     case t_INT:
    1077       72464 :     case t_FRAC: return scalarcol(x, nf_get_degree(nf));
    1078             :   }
    1079           0 :   pari_err_TYPE("algtobasis",x);
    1080             :   return NULL; /* LCOV_EXCL_LINE */
    1081             : }
    1082             : 
    1083             : GEN
    1084       44212 : rnfbasistoalg(GEN rnf,GEN x)
    1085             : {
    1086       44212 :   const char *f = "rnfbasistoalg";
    1087             :   long lx, i;
    1088       44212 :   pari_sp av = avma;
    1089             :   GEN z, nf, R, T;
    1090             : 
    1091       44212 :   checkrnf(rnf);
    1092       44212 :   nf = rnf_get_nf(rnf);
    1093       44212 :   T = nf_get_pol(nf);
    1094       44212 :   R = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
    1095       44212 :   switch(typ(x))
    1096             :   {
    1097             :     case t_COL:
    1098         826 :       z = cgetg_copy(x, &lx);
    1099        2478 :       for (i=1; i<lx; i++)
    1100             :       {
    1101        1701 :         GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
    1102        1652 :         if (typ(c) == t_POL) c = mkpolmod(c,T);
    1103        1652 :         gel(z,i) = c;
    1104             :       }
    1105         777 :       z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
    1106         714 :       return gerepileupto(av, gmodulo(z,R));
    1107             : 
    1108             :     case t_POLMOD:
    1109       29757 :       x = polmod_nffix(f, rnf, x, 0);
    1110       29554 :       if (typ(x) != t_POL) break;
    1111       13545 :       retmkpolmod(RgX_copy(x), RgX_copy(R));
    1112             :     case t_POL:
    1113        1099 :       if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
    1114         875 :       if (varn(x) == varn(R))
    1115             :       {
    1116         826 :         x = RgX_nffix(f,nf_get_pol(nf),x,0);
    1117         826 :         return gmodulo(x, R);
    1118             :       }
    1119          49 :       pari_err_VAR(f, x,R);
    1120             :   }
    1121       28714 :   retmkpolmod(scalarpol(x, varn(R)), RgX_copy(R));
    1122             : }
    1123             : 
    1124             : GEN
    1125        1890 : matbasistoalg(GEN nf,GEN x)
    1126             : {
    1127             :   long i, j, li, lx;
    1128        1890 :   GEN z = cgetg_copy(x, &lx);
    1129             : 
    1130        1890 :   if (lx == 1) return z;
    1131        1883 :   switch(typ(x))
    1132             :   {
    1133             :     case t_VEC: case t_COL:
    1134          56 :       for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
    1135          56 :       return z;
    1136        1827 :     case t_MAT: break;
    1137           0 :     default: pari_err_TYPE("matbasistoalg",x);
    1138             :   }
    1139        1827 :   li = lgcols(x);
    1140        6874 :   for (j=1; j<lx; j++)
    1141             :   {
    1142        5047 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1143        5047 :     gel(z,j) = c;
    1144        5047 :     for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
    1145             :   }
    1146        1827 :   return z;
    1147             : }
    1148             : 
    1149             : GEN
    1150        4634 : matalgtobasis(GEN nf,GEN x)
    1151             : {
    1152             :   long i, j, li, lx;
    1153        4634 :   GEN z = cgetg_copy(x, &lx);
    1154             : 
    1155        4634 :   if (lx == 1) return z;
    1156        4291 :   switch(typ(x))
    1157             :   {
    1158             :     case t_VEC: case t_COL:
    1159        4284 :       for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
    1160        4284 :       return z;
    1161           7 :     case t_MAT: break;
    1162           0 :     default: pari_err_TYPE("matalgtobasis",x);
    1163             :   }
    1164           7 :   li = lgcols(x);
    1165          14 :   for (j=1; j<lx; j++)
    1166             :   {
    1167           7 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1168           7 :     gel(z,j) = c;
    1169           7 :     for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
    1170             :   }
    1171           7 :   return z;
    1172             : }
    1173             : GEN
    1174        9366 : RgM_to_nfM(GEN nf,GEN x)
    1175             : {
    1176             :   long i, j, li, lx;
    1177        9366 :   GEN z = cgetg_copy(x, &lx);
    1178             : 
    1179        9366 :   if (lx == 1) return z;
    1180        9366 :   li = lgcols(x);
    1181       71470 :   for (j=1; j<lx; j++)
    1182             :   {
    1183       62104 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1184       62104 :     gel(z,j) = c;
    1185       62104 :     for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
    1186             :   }
    1187        9366 :   return z;
    1188             : }
    1189             : GEN
    1190       86317 : RgC_to_nfC(GEN nf, GEN x)
    1191       86317 : { pari_APPLY_type(t_COL, nf_to_scalar_or_basis(nf, gel(x,i))) }
    1192             : 
    1193             : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
    1194             : GEN
    1195      134302 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
    1196      134302 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
    1197             : GEN
    1198      134393 : polmod_nffix2(const char *f, GEN T, GEN R, GEN x, int lift)
    1199             : {
    1200      134393 :   if (RgX_equal_var(gel(x,1), R))
    1201             :   {
    1202      124432 :     x = gel(x,2);
    1203      124432 :     if (typ(x) == t_POL && varn(x) == varn(R))
    1204             :     {
    1205       94808 :       x = RgX_nffix(f, T, x, lift);
    1206       94808 :       switch(lg(x))
    1207             :       {
    1208        5782 :         case 2: return gen_0;
    1209       11347 :         case 3: return gel(x,2);
    1210             :       }
    1211       77679 :       return x;
    1212             :     }
    1213             :   }
    1214       39585 :   return Rg_nffix(f, T, x, lift);
    1215             : }
    1216             : GEN
    1217        1204 : rnfalgtobasis(GEN rnf,GEN x)
    1218             : {
    1219        1204 :   const char *f = "rnfalgtobasis";
    1220        1204 :   pari_sp av = avma;
    1221             :   GEN T, R;
    1222             : 
    1223        1204 :   checkrnf(rnf);
    1224        1204 :   R = rnf_get_pol(rnf);
    1225        1204 :   T = rnf_get_nfpol(rnf);
    1226        1204 :   switch(typ(x))
    1227             :   {
    1228             :     case t_COL:
    1229          49 :       if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
    1230          28 :       x = RgV_nffix(f, T, x, 0);
    1231          21 :       return gerepilecopy(av, x);
    1232             : 
    1233             :     case t_POLMOD:
    1234        1071 :       x = polmod_nffix(f, rnf, x, 0);
    1235        1036 :       if (typ(x) != t_POL) break;
    1236         714 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1237             :     case t_POL:
    1238          56 :       if (varn(x) == varn(T))
    1239             :       {
    1240          21 :         RgX_check_QX(x,f);
    1241          14 :         if (degpol(x) >= degpol(T)) x = RgX_rem(x,T);
    1242          14 :         x = mkpolmod(x,T); break;
    1243             :       }
    1244          35 :       x = RgX_nffix(f, T, x, 0);
    1245          28 :       if (degpol(x) >= degpol(R)) x = RgX_rem(x, R);
    1246          28 :       return gerepileupto(av, RgM_RgX_mul(rnf_get_invzk(rnf), x));
    1247             :   }
    1248         364 :   return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
    1249             : }
    1250             : 
    1251             : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
    1252             :  * is "small" */
    1253             : GEN
    1254         259 : nfdiveuc(GEN nf, GEN a, GEN b)
    1255             : {
    1256         259 :   pari_sp av = avma;
    1257         259 :   a = nfdiv(nf,a,b);
    1258         259 :   return gerepileupto(av, ground(a));
    1259             : }
    1260             : 
    1261             : /* Given a and b in nf, gives a "small" algebraic integer r in nf
    1262             :  * of the form a-b.y */
    1263             : GEN
    1264         259 : nfmod(GEN nf, GEN a, GEN b)
    1265             : {
    1266         259 :   pari_sp av = avma;
    1267         259 :   GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
    1268         259 :   return gerepileupto(av, nfadd(nf,a,p1));
    1269             : }
    1270             : 
    1271             : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
    1272             :  * that r=a-b.y is "small". */
    1273             : GEN
    1274         259 : nfdivrem(GEN nf, GEN a, GEN b)
    1275             : {
    1276         259 :   pari_sp av = avma;
    1277         259 :   GEN p1,z, y = ground(nfdiv(nf,a,b));
    1278             : 
    1279         259 :   p1 = gneg_i(nfmul(nf,b,y));
    1280         259 :   z = cgetg(3,t_VEC);
    1281         259 :   gel(z,1) = gcopy(y);
    1282         259 :   gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
    1283             : }
    1284             : 
    1285             : /*************************************************************************/
    1286             : /**                                                                     **/
    1287             : /**                   LOGARITHMIC EMBEDDINGS                            **/
    1288             : /**                                                                     **/
    1289             : /*************************************************************************/
    1290             : 
    1291             : static int
    1292      530233 : low_prec(GEN x)
    1293             : {
    1294      530233 :   switch(typ(x))
    1295             :   {
    1296           0 :     case t_INT: return !signe(x);
    1297      530233 :     case t_REAL: return !signe(x) || realprec(x) <= DEFAULTPREC;
    1298           0 :     default: return 0;
    1299             :   }
    1300             : }
    1301             : 
    1302             : static GEN
    1303       11192 : triv_cxlog(GEN nf) { return zerocol(lg(nf_get_roots(nf))-1); }
    1304             : static GEN
    1305        9583 : famat_cxlog(GEN nf, GEN fa, long prec)
    1306             : {
    1307        9583 :   GEN g,e, y = NULL;
    1308             :   long i,l;
    1309             : 
    1310        9583 :   if (typ(fa) != t_MAT) pari_err_TYPE("famat_cxlog",fa);
    1311        9583 :   if (lg(fa) == 1) return triv_cxlog(nf);
    1312        9583 :   g = gel(fa,1);
    1313        9583 :   e = gel(fa,2); l = lg(e);
    1314       17722 :   for (i = 1; i < l; i++)
    1315             :   {
    1316        8139 :     GEN t, x = nf_to_scalar_or_basis(nf, gel(g,i));
    1317             :     /* multiplicative arch would be better (save logs), but exponents overflow
    1318             :      * [ could keep track of expo separately, but not worth it ] */
    1319        8139 :     t = nf_cxlog(nf,x,prec); if (!t) return NULL;
    1320        8139 :     if (gel(t,1) == gen_0) continue; /* rational */
    1321        4604 :     t = RgC_Rg_mul(t, gel(e,i));
    1322        4604 :     y = y? RgV_add(y,t): t;
    1323             :   }
    1324        9583 :   return y ? y: triv_cxlog(nf);
    1325             : }
    1326             : /* Archimedean components: [e_i Log( sigma_i(X) )], where X = primpart(x),
    1327             :  * and e_i = 1 (resp 2.) for i <= R1 (resp. > R1) */
    1328             : GEN
    1329      229550 : nf_cxlog(GEN nf, GEN x, long prec)
    1330             : {
    1331             :   long i, l, r1;
    1332             :   GEN v;
    1333      229550 :   if (typ(x) == t_MAT) return famat_cxlog(nf,x,prec);
    1334      219967 :   x = nf_to_scalar_or_basis(nf,x);
    1335      219967 :   if (typ(x) != t_COL) return triv_cxlog(nf);
    1336      215788 :   x = RgM_RgC_mul(nf_get_M(nf), Q_primpart(x));
    1337      215788 :   l = lg(x); r1 = nf_get_r1(nf);
    1338      644114 :   for (i = 1; i <= r1; i++)
    1339      428326 :     if (low_prec(gel(x,i))) return NULL;
    1340      309729 :   for (     ; i <  l;  i++)
    1341       93941 :     if (low_prec(gnorm(gel(x,i)))) return NULL;
    1342      215788 :   v = cgetg(l,t_COL);
    1343      215788 :   for (i = 1; i <= r1; i++) gel(v,i) = glog(gel(x,i),prec);
    1344      215788 :   for (     ; i <  l;  i++) gel(v,i) = gmul2n(glog(gel(x,i),prec),1);
    1345      215788 :   return v;
    1346             : }
    1347             : GEN
    1348          91 : nfV_cxlog(GEN nf, GEN x, long prec)
    1349             : {
    1350             :   long i, l;
    1351          91 :   GEN v = cgetg_copy(x, &l);
    1352         161 :   for (i = 1; i < l; i++)
    1353          70 :     if (!(gel(v,i) = nf_cxlog(nf, gel(x,i), prec))) return NULL;
    1354          91 :   return v;
    1355             : }
    1356             : 
    1357             : static GEN
    1358        1477 : scalar_logembed(GEN nf, GEN u, GEN *emb)
    1359             : {
    1360             :   GEN v, logu;
    1361        1477 :   long i, s = signe(u), RU = lg(nf_get_roots(nf))-1, R1 = nf_get_r1(nf);
    1362             : 
    1363        1477 :   if (!s) pari_err_DOMAIN("nflogembed","argument","=",gen_0,u);
    1364        1477 :   v = cgetg(RU+1, t_COL); logu = logr_abs(u);
    1365        1477 :   for (i = 1; i <= R1; i++) gel(v,i) = logu;
    1366        1477 :   if (i <= RU)
    1367             :   {
    1368         609 :     GEN logu2 = shiftr(logu,1);
    1369         609 :     for (   ; i <= RU; i++) gel(v,i) = logu2;
    1370             :   }
    1371        1477 :   if (emb) *emb = const_col(RU, u);
    1372        1477 :   return v;
    1373             : }
    1374             : 
    1375             : static GEN
    1376        1211 : famat_logembed(GEN nf,GEN x,GEN *emb,long prec)
    1377             : {
    1378        1211 :   GEN A, M, T, a, t, g = gel(x,1), e = gel(x,2);
    1379        1211 :   long i, l = lg(e);
    1380             : 
    1381        1211 :   if (l == 1) return scalar_logembed(nf, real_1(prec), emb);
    1382        1211 :   A = NULL; T = emb? cgetg(l, t_COL): NULL;
    1383        1211 :   if (emb) *emb = M = mkmat2(T, e);
    1384        4109 :   for (i = 1; i < l; i++)
    1385             :   {
    1386        2898 :     a = nflogembed(nf, gel(g,i), &t, prec);
    1387        2898 :     if (!a) return NULL;
    1388        2898 :     a = RgC_Rg_mul(a, gel(e,i));
    1389        2898 :     A = A? RgC_add(A, a): a;
    1390        2898 :     if (emb) gel(T,i) = t;
    1391             :   }
    1392        1211 :   return A;
    1393             : }
    1394             : 
    1395             : /* Get archimedean components: [e_i log( | sigma_i(x) | )], with e_i = 1
    1396             :  * (resp 2.) for i <= R1 (resp. > R1) and set emb to the embeddings of x.
    1397             :  * Return NULL if precision problem */
    1398             : GEN
    1399        5712 : nflogembed(GEN nf, GEN x, GEN *emb, long prec)
    1400             : {
    1401             :   long i, l, r1;
    1402             :   GEN v, t;
    1403             : 
    1404        5712 :   if (typ(x) == t_MAT) return famat_logembed(nf,x,emb,prec);
    1405        4501 :   x = nf_to_scalar_or_basis(nf,x);
    1406        4501 :   if (typ(x) != t_COL) return scalar_logembed(nf, gtofp(x,prec), emb);
    1407        3024 :   x = RgM_RgC_mul(nf_get_M(nf), x);
    1408        3024 :   l = lg(x); r1 = nf_get_r1(nf); v = cgetg(l,t_COL);
    1409        7840 :   for (i = 1; i <= r1; i++)
    1410             :   {
    1411        4823 :     t = gabs(gel(x,i),prec); if (low_prec(t)) return NULL;
    1412        4816 :     gel(v,i) = glog(t,prec);
    1413             :   }
    1414        6153 :   for (   ; i < l; i++)
    1415             :   {
    1416        3143 :     t = gnorm(gel(x,i)); if (low_prec(t)) return NULL;
    1417        3136 :     gel(v,i) = glog(t,prec);
    1418             :   }
    1419        3010 :   if (emb) *emb = x;
    1420        3010 :   return v;
    1421             : }
    1422             : 
    1423             : /*************************************************************************/
    1424             : /**                                                                     **/
    1425             : /**                        REAL EMBEDDINGS                              **/
    1426             : /**                                                                     **/
    1427             : /*************************************************************************/
    1428             : static GEN
    1429       83951 : sarch_get_cyc(GEN sarch) { return gel(sarch,1); }
    1430             : static GEN
    1431      340176 : sarch_get_archp(GEN sarch) { return gel(sarch,2); }
    1432             : static GEN
    1433       67007 : sarch_get_MI(GEN sarch) { return gel(sarch,3); }
    1434             : static GEN
    1435       67007 : sarch_get_lambda(GEN sarch) { return gel(sarch,4); }
    1436             : static GEN
    1437       67007 : sarch_get_F(GEN sarch) { return gel(sarch,5); }
    1438             : 
    1439             : /* true nf, return number of positive roots of char_x */
    1440             : static long
    1441        1858 : num_positive(GEN nf, GEN x)
    1442             : {
    1443        1858 :   GEN T = nf_get_pol(nf);
    1444        1858 :   GEN charx = ZXQ_charpoly(nf_to_scalar_or_alg(nf,x), T, 0);
    1445             :   long np;
    1446        1858 :   charx = ZX_radical(charx);
    1447        1858 :   np = ZX_sturmpart(charx, mkvec2(gen_0,mkoo()));
    1448        1858 :   return np * (degpol(T) / degpol(charx));
    1449             : }
    1450             : 
    1451             : /* x a QC: return sigma_k(x) where 1 <= k <= r1+r2; correct but inefficient
    1452             :  * if x in Q. M = nf_get_M(nf) */
    1453             : static GEN
    1454          91 : nfembed_i(GEN M, GEN x, long k)
    1455             : {
    1456          91 :   long i, l = lg(M);
    1457          91 :   GEN z = gel(x,1);
    1458          91 :   for (i = 2; i < l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
    1459          91 :   return z;
    1460             : }
    1461             : GEN
    1462           0 : nfembed(GEN nf, GEN x, long k)
    1463             : {
    1464           0 :   pari_sp av = avma;
    1465           0 :   nf = checknf(nf);
    1466           0 :   x = nf_to_scalar_or_basis(nf,x);
    1467           0 :   if (typ(x) != t_COL) return gerepilecopy(av, x);
    1468           0 :   return gerepileupto(av, nfembed_i(nf_get_M(nf),x,k));
    1469             : }
    1470             : 
    1471             : /* x a ZC */
    1472             : static GEN
    1473      435460 : zk_embed(GEN M, GEN x, long k)
    1474             : {
    1475      435460 :   long i, l = lg(x);
    1476      435460 :   GEN z = gel(x,1); /* times M[k,1], which is 1 */
    1477      435460 :   for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
    1478      435460 :   return z;
    1479             : }
    1480             : 
    1481             : /* Given floating point approximation z of sigma_k(x), decide its sign
    1482             :  * [0/+, 1/- and -1 for FAIL] */
    1483             : static long
    1484      423726 : eval_sign_embed(GEN z)
    1485             : { /* dubious, fail */
    1486      423726 :   if (typ(z) == t_REAL && realprec(z) <= LOWDEFAULTPREC) return -1;
    1487      422487 :   return (signe(z) < 1)? 1: 0;
    1488             : }
    1489             : /* return v such that (-1)^v = sign(sigma_k(x)), x primitive ZC */
    1490             : static long
    1491      346021 : eval_sign(GEN M, GEN x, long k)
    1492      346021 : { return eval_sign_embed( zk_embed(M, x, k) ); }
    1493             : 
    1494             : /* check that signs[i..#signs] == s; signs = NULL encodes "totally positive" */
    1495             : static int
    1496           0 : oksigns(long l, GEN signs, long i, long s)
    1497             : {
    1498           0 :   if (!signs) return s == 0;
    1499           0 :   for (; i < l; i++)
    1500           0 :     if (signs[i] != s) return 0;
    1501           0 :   return 1;
    1502             : }
    1503             : /* check that signs[i] = s and signs[i+1..#signs] = 1-s */
    1504             : static int
    1505           0 : oksigns2(long l, GEN signs, long i, long s)
    1506             : {
    1507           0 :   if (!signs) return s == 0 && i == l-1;
    1508           0 :   return signs[i] == s && oksigns(l, signs, i+1, 1-s);
    1509             : }
    1510             : 
    1511             : /* true nf, x a ZC (primitive for efficiency), embx its embeddings or NULL */
    1512             : static int
    1513       67655 : nfchecksigns_i(GEN nf, GEN x, GEN embx, GEN signs, GEN archp)
    1514             : {
    1515       67655 :   long l = lg(archp), i;
    1516       67655 :   GEN M = nf_get_M(nf), sarch = NULL;
    1517       67655 :   long np = -1;
    1518      102162 :   for (i = 1; i < l; i++)
    1519             :   {
    1520             :     long s;
    1521       78251 :     if (embx)
    1522       77705 :       s = eval_sign_embed(gel(embx,i));
    1523             :     else
    1524         546 :       s = eval_sign(M, x, archp[i]);
    1525             :     /* 0 / + or 1 / -; -1 for FAIL */
    1526       78251 :     if (s < 0) /* failure */
    1527             :     {
    1528           0 :       long ni, r1 = nf_get_r1(nf);
    1529             :       GEN xi;
    1530           0 :       if (np < 0)
    1531             :       {
    1532           0 :         np = num_positive(nf, x);
    1533           0 :         if (np == 0)  return oksigns(l, signs, i, 1);
    1534           0 :         if (np == r1) return oksigns(l, signs, i, 0);
    1535           0 :         sarch = nfarchstar(nf, NULL, identity_perm(r1));
    1536             :       }
    1537           0 :       xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    1538           0 :       xi = Q_primpart(xi);
    1539           0 :       ni = num_positive(nf, nfmuli(nf,x,xi));
    1540           0 :       if (ni == 0)  return oksigns2(l, signs, i, 0);
    1541           0 :       if (ni == r1) return oksigns2(l, signs, i, 1);
    1542           0 :       s = ni < np? 0: 1;
    1543             :     }
    1544       78251 :     if (s != (signs? signs[i]: 0)) return 0;
    1545             :   }
    1546       23911 :   return 1;
    1547             : }
    1548             : static void
    1549         406 : pl_convert(GEN pl, GEN *psigns, GEN *parchp)
    1550             : {
    1551         406 :   long i, j, l = lg(pl);
    1552         406 :   GEN signs = cgetg(l, t_VECSMALL);
    1553         406 :   GEN archp = cgetg(l, t_VECSMALL);
    1554        1372 :   for (i = j = 1; i < l; i++)
    1555             :   {
    1556         966 :     if (!pl[i]) continue;
    1557         714 :     archp[j] = i;
    1558         714 :     signs[j] = (pl[i] < 0)? 1: 0;
    1559         714 :     j++;
    1560             :   }
    1561         406 :   setlg(archp, j); *parchp = archp;
    1562         406 :   setlg(signs, j); *psigns = signs;
    1563         406 : }
    1564             : /* pl : requested signs for real embeddings, 0 = no sign constraint */
    1565             : int
    1566         931 : nfchecksigns(GEN nf, GEN x, GEN pl)
    1567             : {
    1568         931 :   pari_sp av = avma;
    1569             :   GEN signs, archp;
    1570         931 :   nf = checknf(nf);
    1571         931 :   x = nf_to_scalar_or_basis(nf,x);
    1572         931 :   if (typ(x) != t_COL)
    1573             :   {
    1574         525 :     long i, l = lg(pl), s = gsigne(x);
    1575        1064 :     for (i = 1; i < l; i++)
    1576         539 :       if (pl[i] && pl[i] != s) return gc_bool(av,0);
    1577         525 :     return gc_bool(av,1);
    1578             :   }
    1579         406 :   pl_convert(pl, &signs, &archp);
    1580         406 :   return gc_bool(av, nfchecksigns_i(nf, x, NULL, signs, archp));
    1581             : }
    1582             : 
    1583             : /* signs = NULL: totally positive, else sign[i] = 0 (+) or 1 (-) */
    1584             : static GEN
    1585       67007 : get_C(GEN lambda, long l, GEN signs)
    1586             : {
    1587             :   long i;
    1588             :   GEN C, mlambda;
    1589       67007 :   if (!signs) return const_vec(l-1, lambda);
    1590       28325 :   C = cgetg(l, t_COL); mlambda = gneg(lambda);
    1591       28325 :   for (i = 1; i < l; i++) gel(C,i) = signs[i]? mlambda: lambda;
    1592       28325 :   return C;
    1593             : }
    1594             : /* signs = NULL: totally positive at archp */
    1595             : static GEN
    1596      120500 : nfsetsigns(GEN nf, GEN signs, GEN x, GEN sarch)
    1597             : {
    1598      120500 :   long i, l = lg(sarch_get_archp(sarch));
    1599             :   GEN ex;
    1600             :   /* Is signature already correct ? */
    1601      120500 :   if (typ(x) != t_COL)
    1602             :   {
    1603       53251 :     long s = gsigne(x);
    1604       53251 :     if (!s) i = 1;
    1605       53237 :     else if (!signs)
    1606        4942 :       i = (s < 0)? 1: l;
    1607             :     else
    1608             :     {
    1609       48295 :       s = s < 0? 1: 0;
    1610       78207 :       for (i = 1; i < l; i++)
    1611       53154 :         if (signs[i] != s) break;
    1612             :     }
    1613       53251 :     ex = (i < l)? const_col(l-1, x): NULL;
    1614             :   }
    1615             :   else
    1616             :   {
    1617       67249 :     pari_sp av = avma;
    1618       67249 :     GEN cex, M = nf_get_M(nf), archp = sarch_get_archp(sarch);
    1619       67249 :     GEN xp = Q_primitive_part(x,&cex);
    1620       67249 :     ex = cgetg(l,t_COL);
    1621       67249 :     for (i = 1; i < l; i++) gel(ex,i) = zk_embed(M,xp,archp[i]);
    1622       67249 :     if (nfchecksigns_i(nf, xp, ex, signs, archp)) { ex = NULL; set_avma(av); }
    1623       43625 :     else if (cex) ex = RgC_Rg_mul(ex, cex); /* put back content */
    1624             :   }
    1625      120500 :   if (ex)
    1626             :   { /* If no, fix it */
    1627       67007 :     GEN MI = sarch_get_MI(sarch), F = sarch_get_F(sarch);
    1628       67007 :     GEN lambda = sarch_get_lambda(sarch);
    1629       67007 :     GEN t = RgC_sub(get_C(lambda, l, signs), ex);
    1630             :     long e;
    1631       67007 :     t = grndtoi(RgM_RgC_mul(MI,t), &e);
    1632       67007 :     if (lg(F) != 1) t = ZM_ZC_mul(F, t);
    1633       67007 :     x = typ(x) == t_COL? RgC_add(t, x): RgC_Rg_add(t, x);
    1634             :   }
    1635      120500 :   return x;
    1636             : }
    1637             : /* - sarch = nfarchstar(nf, F);
    1638             :  * - x encodes a vector of signs at arch.archp: either a t_VECSMALL
    1639             :  *   (vector of signs as {0,1}-vector), NULL (totally positive at archp),
    1640             :  *   or a non-zero number field element (replaced by its signature at archp);
    1641             :  * - y is a non-zero number field element
    1642             :  * Return z = y (mod F) with signs(y, archp) = signs(x) (a {0,1}-vector) */
    1643             : GEN
    1644      152427 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN sarch)
    1645             : {
    1646      152427 :   GEN archp = sarch_get_archp(sarch);
    1647      152427 :   if (lg(archp) == 1) return y;
    1648      119135 :   nf = checknf(nf);
    1649      119135 :   if (x && typ(x) != t_VECSMALL) x = nfsign_arch(nf, x, archp);
    1650      119135 :   y = nf_to_scalar_or_basis(nf,y);
    1651      119135 :   return nfsetsigns(nf, x, y, sarch);
    1652             : }
    1653             : 
    1654             : static GEN
    1655       23642 : setsigns_init(GEN nf, GEN archp, GEN F, GEN DATA)
    1656             : {
    1657       23642 :   GEN lambda, Mr = rowpermute(nf_get_M(nf), archp), MI = F? RgM_mul(Mr,F): Mr;
    1658       23642 :   lambda = gmul2n(matrixnorm(MI,DEFAULTPREC), -1);
    1659       23642 :   if (typ(lambda) != t_REAL) lambda = gmul(lambda, sstoQ(1001,1000));
    1660       23642 :   if (lg(archp) < lg(MI))
    1661             :   {
    1662       21203 :     GEN perm = gel(indexrank(MI), 2);
    1663       21203 :     if (!F) F = matid(nf_get_degree(nf));
    1664       21203 :     MI = vecpermute(MI, perm);
    1665       21203 :     F = vecpermute(F, perm);
    1666             :   }
    1667       23642 :   if (!F) F = cgetg(1,t_MAT);
    1668       23642 :   MI = RgM_inv(MI);
    1669       23642 :   return mkvec5(DATA, archp, MI, lambda, F);
    1670             : }
    1671             : /* F non-0 integral ideal in HNF (or NULL: Z_K), compute elements in 1+F
    1672             :  * whose sign matrix at archp is identity; archp in 'indices' format */
    1673             : GEN
    1674       46931 : nfarchstar(GEN nf, GEN F, GEN archp)
    1675             : {
    1676       46931 :   long nba = lg(archp) - 1;
    1677       46931 :   if (!nba) return mkvec2(cgetg(1,t_VEC), archp);
    1678       22284 :   if (F && equali1(gcoeff(F,1,1))) F = NULL;
    1679       22284 :   if (F) F = idealpseudored(F, nf_get_roundG(nf));
    1680       22284 :   return setsigns_init(nf, archp, F, const_vec(nba, gen_2));
    1681             : }
    1682             : 
    1683             : /*************************************************************************/
    1684             : /**                                                                     **/
    1685             : /**                         IDEALCHINESE                                **/
    1686             : /**                                                                     **/
    1687             : /*************************************************************************/
    1688             : static int
    1689        2982 : isprfact(GEN x)
    1690             : {
    1691             :   long i, l;
    1692             :   GEN L, E;
    1693        2982 :   if (typ(x) != t_MAT || lg(x) != 3) return 0;
    1694        2982 :   L = gel(x,1); l = lg(L);
    1695        2982 :   E = gel(x,2);
    1696        7133 :   for(i=1; i<l; i++)
    1697             :   {
    1698        4151 :     checkprid(gel(L,i));
    1699        4151 :     if (typ(gel(E,i)) != t_INT) return 0;
    1700             :   }
    1701        2982 :   return 1;
    1702             : }
    1703             : 
    1704             : /* initialize projectors mod pr[i]^e[i] for idealchinese */
    1705             : static GEN
    1706        2982 : pr_init(GEN nf, GEN fa, GEN w, GEN dw)
    1707             : {
    1708        2982 :   GEN U, E, F, L = gel(fa,1), E0 = gel(fa,2);
    1709        2982 :   long i, r = lg(L);
    1710             : 
    1711        2982 :   if (w && lg(w) != r) pari_err_TYPE("idealchinese", w);
    1712        2982 :   if (r == 1 && !dw) return cgetg(1,t_VEC);
    1713        2975 :   E = leafcopy(E0); /* do not destroy fa[2] */
    1714        7126 :   for (i = 1; i < r; i++)
    1715        4151 :     if (signe(gel(E,i)) < 0) gel(E,i) = gen_0;
    1716        2975 :   F = factorbackprime(nf, L, E);
    1717        2975 :   if (dw)
    1718             :   {
    1719         700 :     F = ZM_Z_mul(F, dw);
    1720        1603 :     for (i = 1; i < r; i++)
    1721             :     {
    1722         903 :       GEN pr = gel(L,i);
    1723         903 :       long e = itos(gel(E0,i)), v = idealval(nf, dw, pr);
    1724         903 :       if (e >= 0)
    1725         896 :         gel(E,i) = addiu(gel(E,i), v);
    1726           7 :       else if (v + e <= 0)
    1727           0 :         F = idealmulpowprime(nf, F, pr, stoi(-v)); /* coprime to pr */
    1728             :       else
    1729             :       {
    1730           7 :         F = idealmulpowprime(nf, F, pr, stoi(e));
    1731           7 :         gel(E,i) = stoi(v + e);
    1732             :       }
    1733             :     }
    1734             :   }
    1735        2975 :   U = cgetg(r, t_VEC);
    1736        7126 :   for (i = 1; i < r; i++)
    1737             :   {
    1738             :     GEN u;
    1739        4151 :     if (w && gequal0(gel(w,i))) u = gen_0; /* unused */
    1740             :     else
    1741             :     {
    1742        4074 :       GEN pr = gel(L,i), e = gel(E,i), t;
    1743        4074 :       t = idealdivpowprime(nf,F, pr, e);
    1744        4074 :       u = hnfmerge_get_1(t, idealpow(nf, pr, e));
    1745        4074 :       if (!u) pari_err_COPRIME("idealchinese", t,pr);
    1746             :     }
    1747        4151 :     gel(U,i) = u;
    1748             :   }
    1749        2975 :   F = idealpseudored(F, nf_get_roundG(nf));
    1750        2975 :   return mkvec2(F, U);
    1751             : }
    1752             : 
    1753             : static GEN
    1754        1778 : pl_normalize(GEN nf, GEN pl)
    1755             : {
    1756        1778 :   const char *fun = "idealchinese";
    1757        1778 :   if (lg(pl)-1 != nf_get_r1(nf)) pari_err_TYPE(fun,pl);
    1758        1778 :   switch(typ(pl))
    1759             :   {
    1760         714 :     case t_VEC: RgV_check_ZV(pl,fun); pl = ZV_to_zv(pl);
    1761             :       /* fall through */
    1762        1778 :     case t_VECSMALL: break;
    1763           0 :     default: pari_err_TYPE(fun,pl);
    1764             :   }
    1765        1778 :   return pl;
    1766             : }
    1767             : 
    1768             : static int
    1769        7070 : is_chineseinit(GEN x)
    1770             : {
    1771             :   GEN fa, pl;
    1772             :   long l;
    1773        7070 :   if (typ(x) != t_VEC || lg(x)!=3) return 0;
    1774        5397 :   fa = gel(x,1);
    1775        5397 :   pl = gel(x,2);
    1776        5397 :   if (typ(fa) != t_VEC || typ(pl) != t_VEC) return 0;
    1777        2639 :   l = lg(fa);
    1778        2639 :   if (l != 1)
    1779             :   {
    1780        2618 :     if (l != 3 || typ(gel(fa,1)) != t_MAT || typ(gel(fa,2)) != t_VEC)
    1781           7 :       return 0;
    1782             :   }
    1783        2632 :   l = lg(pl);
    1784        2632 :   if (l != 1)
    1785             :   {
    1786         511 :     if (l != 6 || typ(gel(pl,3)) != t_MAT || typ(gel(pl,1)) != t_VECSMALL
    1787         511 :                                           || typ(gel(pl,2)) != t_VECSMALL)
    1788           0 :       return 0;
    1789             :   }
    1790        2632 :   return 1;
    1791             : }
    1792             : 
    1793             : /* nf a true 'nf' */
    1794             : static GEN
    1795        3115 : chineseinit_i(GEN nf, GEN fa, GEN w, GEN dw)
    1796             : {
    1797        3115 :   const char *fun = "idealchineseinit";
    1798        3115 :   GEN archp = NULL, pl = NULL;
    1799        3115 :   switch(typ(fa))
    1800             :   {
    1801             :     case t_VEC:
    1802        1778 :       if (is_chineseinit(fa))
    1803             :       {
    1804           0 :         if (dw) pari_err_DOMAIN(fun, "denom(y)", "!=", gen_1, w);
    1805           0 :         return fa;
    1806             :       }
    1807        1778 :       if (lg(fa) != 3) pari_err_TYPE(fun, fa);
    1808             :       /* of the form [x,s] */
    1809        1778 :       pl = pl_normalize(nf, gel(fa,2));
    1810        1778 :       fa = gel(fa,1);
    1811        1778 :       archp = vecsmall01_to_indices(pl);
    1812             :       /* keep pr_init, reset pl */
    1813        1778 :       if (is_chineseinit(fa)) { fa = gel(fa,1); break; }
    1814             :       /* fall through */
    1815             :     case t_MAT: /* factorization? */
    1816        2982 :       if (isprfact(fa)) { fa = pr_init(nf, fa, w, dw); break; }
    1817           0 :     default: pari_err_TYPE(fun,fa);
    1818             :   }
    1819             : 
    1820        3115 :   if (!pl) pl = cgetg(1,t_VEC);
    1821             :   else
    1822             :   {
    1823        1778 :     long r = lg(archp);
    1824        1778 :     if (r == 1) pl = cgetg(1, t_VEC);
    1825             :     else
    1826             :     {
    1827        1358 :       GEN F = (lg(fa) == 1)? NULL: gel(fa,1), signs = cgetg(r, t_VECSMALL);
    1828             :       long i;
    1829        1358 :       for (i = 1; i < r; i++) signs[i] = (pl[archp[i]] < 0)? 1: 0;
    1830        1358 :       pl = setsigns_init(nf, archp, F, signs);
    1831             :     }
    1832             :   }
    1833        3115 :   return mkvec2(fa, pl);
    1834             : }
    1835             : 
    1836             : /* Given a prime ideal factorization x, possibly with 0 or negative exponents,
    1837             :  * and a vector w of elements of nf, gives b such that
    1838             :  * v_p(b-w_p)>=v_p(x) for all prime ideals p in the ideal factorization
    1839             :  * and v_p(b)>=0 for all other p, using the standard proof given in GTM 138. */
    1840             : GEN
    1841        5614 : idealchinese(GEN nf, GEN x, GEN w)
    1842             : {
    1843        5614 :   const char *fun = "idealchinese";
    1844        5614 :   pari_sp av = avma;
    1845             :   GEN x1, x2, s, dw, F;
    1846             : 
    1847        5614 :   nf = checknf(nf);
    1848        5614 :   if (!w) return gerepilecopy(av, chineseinit_i(nf,x,NULL,NULL));
    1849             : 
    1850        3514 :   if (typ(w) != t_VEC) pari_err_TYPE(fun,w);
    1851        3514 :   w = Q_remove_denom(matalgtobasis(nf,w), &dw);
    1852        3514 :   if (!is_chineseinit(x)) x = chineseinit_i(nf,x,w,dw);
    1853             :   /* x is a 'chineseinit' */
    1854        3514 :   x1 = gel(x,1); s = NULL;
    1855        3514 :   x2 = gel(x,2);
    1856        3514 :   if (lg(x1) == 1) F = NULL;
    1857             :   else
    1858             :   {
    1859        3493 :     GEN  U = gel(x1,2);
    1860        3493 :     long i, r = lg(w);
    1861        3493 :     F = gel(x1,1);
    1862        9954 :     for (i=1; i<r; i++)
    1863        6461 :       if (!gequal0(gel(w,i)))
    1864             :       {
    1865        4095 :         GEN t = nfmuli(nf, gel(U,i), gel(w,i));
    1866        4095 :         s = s? ZC_add(s,t): t;
    1867             :       }
    1868        3493 :     if (s) s = ZC_reducemodmatrix(s, F);
    1869             :   }
    1870        3514 :   if (lg(x2) != 1) s = nfsetsigns(nf, gel(x2,1), s? s: gen_0, x2);
    1871        3514 :   if (!s) { s = zerocol(nf_get_degree(nf)); dw = NULL; }
    1872             : 
    1873        3514 :   if (dw) s = RgC_Rg_div(s,dw);
    1874        3514 :   return gerepileupto(av, s);
    1875             : }
    1876             : 
    1877             : /*************************************************************************/
    1878             : /**                                                                     **/
    1879             : /**                           (Z_K/I)^*                                 **/
    1880             : /**                                                                     **/
    1881             : /*************************************************************************/
    1882             : GEN
    1883        1778 : vecsmall01_to_indices(GEN v)
    1884             : {
    1885        1778 :   long i, k, l = lg(v);
    1886        1778 :   GEN p = new_chunk(l) + l;
    1887        4774 :   for (k=1, i=l-1; i; i--)
    1888        2996 :     if (v[i]) { *--p = i; k++; }
    1889        1778 :   *--p = evallg(k) | evaltyp(t_VECSMALL);
    1890        1778 :   set_avma((pari_sp)p); return p;
    1891             : }
    1892             : GEN
    1893      404780 : vec01_to_indices(GEN v)
    1894             : {
    1895             :   long i, k, l;
    1896             :   GEN p;
    1897             : 
    1898      404780 :   switch (typ(v))
    1899             :   {
    1900      380980 :    case t_VECSMALL: return v;
    1901       23800 :    case t_VEC: break;
    1902           0 :    default: pari_err_TYPE("vec01_to_indices",v);
    1903             :   }
    1904       23800 :   l = lg(v);
    1905       23800 :   p = new_chunk(l) + l;
    1906       70378 :   for (k=1, i=l-1; i; i--)
    1907       46578 :     if (signe(gel(v,i))) { *--p = i; k++; }
    1908       23800 :   *--p = evallg(k) | evaltyp(t_VECSMALL);
    1909       23800 :   set_avma((pari_sp)p); return p;
    1910             : }
    1911             : GEN
    1912        6545 : indices_to_vec01(GEN p, long r)
    1913             : {
    1914        6545 :   long i, l = lg(p);
    1915        6545 :   GEN v = zerovec(r);
    1916        6545 :   for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
    1917        6545 :   return v;
    1918             : }
    1919             : 
    1920             : /* return (column) vector of R1 signatures of x (0 or 1) */
    1921             : GEN
    1922      380980 : nfsign_arch(GEN nf, GEN x, GEN arch)
    1923             : {
    1924      380980 :   GEN sarch, M, V, archp = vec01_to_indices(arch);
    1925      380980 :   long i, s, np, n = lg(archp)-1;
    1926             :   pari_sp av;
    1927             : 
    1928      380980 :   if (!n) return cgetg(1,t_VECSMALL);
    1929      378824 :   nf = checknf(nf);
    1930      378824 :   if (typ(x) == t_MAT)
    1931             :   { /* factorisation */
    1932      105949 :     GEN g = gel(x,1), e = gel(x,2);
    1933      105949 :     long l = lg(g);
    1934      105949 :     V = zero_zv(n);
    1935      311463 :     for (i = 1; i < l; i++)
    1936      205514 :       if (mpodd(gel(e,i)))
    1937      186571 :         Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
    1938      105949 :     set_avma((pari_sp)V); return V;
    1939             :   }
    1940      272875 :   av = avma; V = cgetg(n+1,t_VECSMALL);
    1941      272875 :   x = nf_to_scalar_or_basis(nf, x);
    1942      272875 :   switch(typ(x))
    1943             :   {
    1944             :     case t_INT:
    1945       89177 :       s = signe(x);
    1946       89177 :       if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
    1947       89177 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    1948             :     case t_FRAC:
    1949          84 :       s = signe(gel(x,1));
    1950          84 :       set_avma(av); return const_vecsmall(n, (s < 0)? 1: 0);
    1951             :   }
    1952      183614 :   x = Q_primpart(x); M = nf_get_M(nf); sarch = NULL; np = -1;
    1953      527891 :   for (i = 1; i <= n; i++)
    1954             :   {
    1955      345475 :     long s = eval_sign(M, x, archp[i]);
    1956      345475 :     if (s < 0) /* failure */
    1957             :     {
    1958        1239 :       long ni, r1 = nf_get_r1(nf);
    1959             :       GEN xi;
    1960        1239 :       if (np < 0)
    1961             :       {
    1962        1219 :         np = num_positive(nf, x);
    1963        1219 :         if (np == 0) { set_avma(av); return const_vecsmall(n, 1); }
    1964        1074 :         if (np == r1){ set_avma(av); return const_vecsmall(n, 0); }
    1965         619 :         sarch = nfarchstar(nf, NULL, identity_perm(r1));
    1966             :       }
    1967         639 :       xi = set_sign_mod_divisor(nf, vecsmall_ei(r1, archp[i]), gen_1, sarch);
    1968         639 :       xi = Q_primpart(xi);
    1969         639 :       ni = num_positive(nf, nfmuli(nf,x,xi));
    1970         639 :       if (ni == 0) { set_avma(av); V = const_vecsmall(n, 1); V[i] = 0; return V; }
    1971         447 :       if (ni == r1){ set_avma(av); V = const_vecsmall(n, 0); V[i] = 1; return V; }
    1972          41 :       s = ni < np? 0: 1;
    1973             :     }
    1974      344277 :     V[i] = s;
    1975             :   }
    1976      182416 :   set_avma((pari_sp)V); return V;
    1977             : }
    1978             : static void
    1979        6909 : chk_ind(const char *s, long i, long r1)
    1980             : {
    1981        6909 :   if (i <= 0) pari_err_DOMAIN(s, "index", "<=", gen_0, stoi(i));
    1982        6895 :   if (i > r1) pari_err_DOMAIN(s, "index", ">", utoi(r1), utoi(i));
    1983        6860 : }
    1984             : static GEN
    1985        6363 : parse_embed(GEN ind, long r, const char *f)
    1986             : {
    1987             :   long l, i;
    1988        6363 :   if (!ind) return identity_perm(r);
    1989        4893 :   switch(typ(ind))
    1990             :   {
    1991         875 :     case t_INT: case t_VEC: case t_COL: ind = gtovecsmall(ind); break;
    1992        4018 :     case t_VECSMALL: break;
    1993           0 :     default: pari_err_TYPE(f, ind);
    1994             :   }
    1995        4893 :   l = lg(ind);
    1996        4893 :   for (i = 1; i < l; i++) chk_ind(f, ind[i], r);
    1997        4844 :   return ind;
    1998             : }
    1999             : GEN
    2000        4809 : nfeltsign(GEN nf, GEN x, GEN ind0)
    2001             : {
    2002        4809 :   pari_sp av = avma;
    2003             :   long i, l;
    2004             :   GEN v, ind;
    2005        4809 :   nf = checknf(nf);
    2006        4809 :   ind = parse_embed(ind0, nf_get_r1(nf), "nfeltsign");
    2007        4788 :   l = lg(ind);
    2008        4788 :   if (is_rational_t(typ(x)))
    2009             :   { /* nfsign_arch would test this, but avoid converting t_VECSMALL -> t_VEC */
    2010             :     GEN s;
    2011        2156 :     switch(gsigne(x))
    2012             :     {
    2013         532 :       case -1:s = gen_m1; break;
    2014        1617 :       case 1: s = gen_1; break;
    2015           7 :       default: s = gen_0; break;
    2016             :     }
    2017        2156 :     set_avma(av);
    2018        2156 :     return (ind0 && typ(ind0) == t_INT)? s: const_vec(l-1, s);
    2019             :   }
    2020        2632 :   v = nfsign_arch(nf, x, ind);
    2021        2632 :   if (ind0 && typ(ind0) == t_INT) { set_avma(av); return v[1]? gen_m1: gen_1; }
    2022        2618 :   settyp(v, t_VEC);
    2023        2618 :   for (i = 1; i < l; i++) gel(v,i) = v[i]? gen_m1: gen_1;
    2024        2618 :   return gerepileupto(av, v);
    2025             : }
    2026             : 
    2027             : GEN
    2028          63 : nfeltembed(GEN nf, GEN x, GEN ind0, long prec0)
    2029             : {
    2030          63 :   pari_sp av = avma;
    2031             :   long i, e, l, r1, r2, prec, prec1;
    2032             :   GEN v, ind, cx;
    2033          63 :   nf = checknf(nf); nf_get_sign(nf,&r1,&r2);
    2034          63 :   x = nf_to_scalar_or_basis(nf, x);
    2035          56 :   ind = parse_embed(ind0, r1+r2, "nfeltembed");
    2036          49 :   l = lg(ind);
    2037          49 :   if (typ(x) != t_COL)
    2038             :   {
    2039           0 :     if (!(ind0 && typ(ind0) == t_INT)) x = const_vec(l-1, x);
    2040           0 :     return gerepilecopy(av, x);
    2041             :   }
    2042          49 :   x = Q_primitive_part(x, &cx);
    2043          49 :   prec1 = prec0; e = gexpo(x);
    2044          49 :   if (e > 8) prec1 += nbits2extraprec(e);
    2045          49 :   prec = prec1;
    2046          49 :   if (nf_get_prec(nf) < prec) nf = nfnewprec_shallow(nf, prec);
    2047          49 :   v = cgetg(l, t_VEC);
    2048             :   for(;;)
    2049           7 :   {
    2050          56 :     GEN M = nf_get_M(nf);
    2051         140 :     for (i = 1; i < l; i++)
    2052             :     {
    2053          91 :       GEN t = nfembed_i(M, x, ind[i]);
    2054          91 :       long e = gexpo(t);
    2055          91 :       if (gequal0(t) || precision(t) < prec0
    2056          91 :                      || (e < 0 && prec < prec1 + nbits2extraprec(-e)) ) break;
    2057          84 :       if (cx) t = gmul(t, cx);
    2058          84 :       gel(v,i) = t;
    2059             :     }
    2060          56 :     if (i == l) break;
    2061           7 :     prec = precdbl(prec);
    2062           7 :     if (DEBUGLEVEL>1) pari_warn(warnprec,"eltnfembed", prec);
    2063           7 :     nf = nfnewprec_shallow(nf, prec);
    2064             :   }
    2065          49 :   if (ind0 && typ(ind0) == t_INT) v = gel(v,1);
    2066          49 :   return gerepilecopy(av, v);
    2067             : }
    2068             : 
    2069             : /* number of distinct roots of sigma(f) */
    2070             : GEN
    2071        1498 : nfpolsturm(GEN nf, GEN f, GEN ind0)
    2072             : {
    2073        1498 :   pari_sp av = avma;
    2074             :   long d, l, r1, single;
    2075             :   GEN ind, u, v, vr1, T, s, t;
    2076             : 
    2077        1498 :   nf = checknf(nf); T = nf_get_pol(nf); r1 = nf_get_r1(nf);
    2078        1498 :   ind = parse_embed(ind0, r1, "nfpolsturm");
    2079        1477 :   single = ind0 && typ(ind0) == t_INT;
    2080        1477 :   l = lg(ind);
    2081             : 
    2082        1477 :   if (gequal0(f)) pari_err_ROOTS0("nfpolsturm");
    2083        1470 :   if (typ(f) == t_POL && varn(f) != varn(T))
    2084             :   {
    2085        1449 :     f = RgX_nffix("nfsturn", T, f,1);
    2086        1449 :     if (lg(f) == 3) f = NULL;
    2087             :   }
    2088             :   else
    2089             :   {
    2090          21 :     (void)Rg_nffix("nfpolsturm", T, f, 0);
    2091          21 :     f = NULL;
    2092             :   }
    2093        1470 :   if (!f) { set_avma(av); return single? gen_0: zerovec(l-1); }
    2094        1449 :   d = degpol(f);
    2095        1449 :   if (d == 1) { set_avma(av); return single? gen_1: const_vec(l-1,gen_1); }
    2096             : 
    2097        1428 :   vr1 = const_vecsmall(l-1, 1);
    2098        1428 :   u = Q_primpart(f); s = ZV_to_zv(nfeltsign(nf, gel(u,d+2), ind));
    2099        1428 :   v = RgX_deriv(u); t = odd(d)? leafcopy(s): zv_neg(s);
    2100             :   for(;;)
    2101         154 :   {
    2102        1582 :     GEN r = RgX_neg( Q_primpart(RgX_pseudorem(u, v)) ), sr;
    2103        1582 :     long i, dr = degpol(r);
    2104        1582 :     if (dr < 0) break;
    2105        1582 :     sr = ZV_to_zv(nfeltsign(nf, gel(r,dr+2), ind));
    2106        3941 :     for (i = 1; i < l; i++)
    2107        2359 :       if (sr[i] != s[i]) { s[i] = sr[i], vr1[i]--; }
    2108        1582 :     if (odd(dr)) sr = zv_neg(sr);
    2109        3941 :     for (i = 1; i < l; i++)
    2110        2359 :       if (sr[i] != t[i]) { t[i] = sr[i], vr1[i]++; }
    2111        1582 :     if (!dr) break;
    2112         154 :     u = v; v = r;
    2113             :   }
    2114        1428 :   if (single) { set_avma(av); return stoi(vr1[1]); }
    2115         721 :   return gerepileupto(av, zv_to_ZV(vr1));
    2116             : }
    2117             : 
    2118             : 
    2119             : /* return the vector of signs of x; the matrix of such if x is a vector
    2120             :  * of nf elements */
    2121             : GEN
    2122        2198 : nfsign(GEN nf, GEN x)
    2123             : {
    2124             :   long i, l;
    2125             :   GEN archp, S;
    2126             : 
    2127        2198 :   nf = checknf(nf);
    2128        2198 :   archp = identity_perm( nf_get_r1(nf) );
    2129        2198 :   if (typ(x) != t_VEC) return nfsign_arch(nf, x, archp);
    2130         994 :   l = lg(x); S = cgetg(l, t_MAT);
    2131         994 :   for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), archp);
    2132         994 :   return S;
    2133             : }
    2134             : 
    2135             : /* x integral elt, A integral ideal in HNF; reduce x mod A */
    2136             : static GEN
    2137     1748758 : zk_modHNF(GEN x, GEN A)
    2138     1748758 : { return (typ(x) == t_COL)?  ZC_hnfrem(x, A): modii(x, gcoeff(A,1,1)); }
    2139             : 
    2140             : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
    2141             :    outputs an element inverse of x modulo y */
    2142             : GEN
    2143         441 : nfinvmodideal(GEN nf, GEN x, GEN y)
    2144             : {
    2145         441 :   pari_sp av = avma;
    2146         441 :   GEN a, yZ = gcoeff(y,1,1);
    2147             : 
    2148         441 :   if (equali1(yZ)) return gen_0;
    2149         441 :   x = nf_to_scalar_or_basis(nf, x);
    2150         441 :   if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
    2151             : 
    2152         371 :   a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
    2153         371 :   if (!a) pari_err_INV("nfinvmodideal", x);
    2154         371 :   return gerepileupto(av, zk_modHNF(nfdiv(nf,a,x), y));
    2155             : }
    2156             : 
    2157             : static GEN
    2158      912763 : nfsqrmodideal(GEN nf, GEN x, GEN id)
    2159      912763 : { return zk_modHNF(nfsqri(nf,x), id); }
    2160             : static GEN
    2161     1440496 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id)
    2162     1440496 : { return x? zk_modHNF(nfmuli(nf,x,y), id): y; }
    2163             : /* assume x integral, k integer, A in HNF */
    2164             : GEN
    2165      755557 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
    2166             : {
    2167      755557 :   long s = signe(k);
    2168             :   pari_sp av;
    2169             :   GEN y;
    2170             : 
    2171      755557 :   if (!s) return gen_1;
    2172      755557 :   av = avma;
    2173      755557 :   x = nf_to_scalar_or_basis(nf, x);
    2174      755557 :   if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
    2175      368227 :   if (s < 0) { x = nfinvmodideal(nf, x,A); k = negi(k); }
    2176      368227 :   for(y = NULL;;)
    2177             :   {
    2178     2193753 :     if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
    2179     1280990 :     k = shifti(k,-1); if (!signe(k)) break;
    2180      912763 :     x = nfsqrmodideal(nf,x,A);
    2181             :   }
    2182      368227 :   return gerepileupto(av, y);
    2183             : }
    2184             : 
    2185             : /* a * g^n mod id */
    2186             : static GEN
    2187      689701 : nfmulpowmodideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
    2188             : {
    2189      689701 :   return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
    2190             : }
    2191             : 
    2192             : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
    2193             :  * EX = multiple of exponent of (O_K/id)^* */
    2194             : GEN
    2195      259647 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
    2196             : {
    2197      259647 :   GEN EXo2, plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
    2198      259647 :   long i, lx = lg(g);
    2199             : 
    2200      259647 :   if (equali1(idZ)) return gen_1; /* id = Z_K */
    2201      259192 :   EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
    2202     1227993 :   for (i = 1; i < lx; i++)
    2203             :   {
    2204      968801 :     GEN h, n = centermodii(gel(e,i), EX, EXo2);
    2205      968801 :     long sn = signe(n);
    2206      968801 :     if (!sn) continue;
    2207             : 
    2208      587736 :     h = nf_to_scalar_or_basis(nf, gel(g,i));
    2209      587736 :     switch(typ(h))
    2210             :     {
    2211      330290 :       case t_INT: break;
    2212             :       case t_FRAC:
    2213           0 :         h = Fp_div(gel(h,1), gel(h,2), idZ); break;
    2214             :       default:
    2215             :       {
    2216             :         GEN dh;
    2217      257446 :         h = Q_remove_denom(h, &dh);
    2218      257446 :         if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
    2219             :       }
    2220             :     }
    2221      587736 :     if (sn > 0)
    2222      557174 :       plus = nfmulpowmodideal(nf, plus, h, n, id);
    2223             :     else /* sn < 0 */
    2224       30562 :       minus = nfmulpowmodideal(nf, minus, h, negi(n), id);
    2225             :   }
    2226      259192 :   if (minus) plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
    2227      259192 :   return plus? plus: gen_1;
    2228             : }
    2229             : 
    2230             : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute (1+x)/(1+y) in
    2231             :  * the form [[cyc],[gen], U], where U := ux^-1 as a pair [ZM, denom(U)] */
    2232             : static GEN
    2233       29897 : zidealij(GEN x, GEN y)
    2234             : {
    2235       29897 :   GEN U, G, cyc, xp = gcoeff(x,1,1), xi = hnf_invscale(x, xp);
    2236             :   long j, N;
    2237             : 
    2238             :   /* x^(-1) y = relations between the 1 + x_i (HNF) */
    2239       29897 :   cyc = ZM_snf_group(ZM_Z_divexact(ZM_mul(xi, y), xp), &U, &G);
    2240       29897 :   N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
    2241       98427 :   for (j=1; j<N; j++)
    2242             :   {
    2243       68530 :     GEN c = gel(G,j);
    2244       68530 :     gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
    2245       68530 :     if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
    2246             :   }
    2247       29897 :   return mkvec4(cyc, G, ZM_mul(U,xi), xp);
    2248             : }
    2249             : 
    2250             : /* lg(x) > 1, x + 1; shallow */
    2251             : static GEN
    2252       12929 : ZC_add1(GEN x)
    2253             : {
    2254       12929 :   long i, l = lg(x);
    2255       12929 :   GEN y = cgetg(l, t_COL);
    2256       12929 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2257       12929 :   gel(y,1) = addiu(gel(x,1), 1); return y;
    2258             : }
    2259             : /* lg(x) > 1, x - 1; shallow */
    2260             : static GEN
    2261        6363 : ZC_sub1(GEN x)
    2262             : {
    2263        6363 :   long i, l = lg(x);
    2264        6363 :   GEN y = cgetg(l, t_COL);
    2265        6363 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    2266        6363 :   gel(y,1) = subiu(gel(x,1), 1); return y;
    2267             : }
    2268             : 
    2269             : /* x,y are t_INT or ZC */
    2270             : static GEN
    2271           0 : zkadd(GEN x, GEN y)
    2272             : {
    2273           0 :   long tx = typ(x);
    2274           0 :   if (tx == typ(y))
    2275           0 :     return tx == t_INT? addii(x,y): ZC_add(x,y);
    2276             :   else
    2277           0 :     return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
    2278             : }
    2279             : /* x a t_INT or ZC, x+1; shallow */
    2280             : static GEN
    2281       27321 : zkadd1(GEN x)
    2282             : {
    2283       27321 :   long tx = typ(x);
    2284       27321 :   return tx == t_INT? addiu(x,1): ZC_add1(x);
    2285             : }
    2286             : /* x a t_INT or ZC, x-1; shallow */
    2287             : static GEN
    2288       27321 : zksub1(GEN x)
    2289             : {
    2290       27321 :   long tx = typ(x);
    2291       27321 :   return tx == t_INT? subiu(x,1): ZC_sub1(x);
    2292             : }
    2293             : /* x,y are t_INT or ZC; x - y */
    2294             : static GEN
    2295           0 : zksub(GEN x, GEN y)
    2296             : {
    2297           0 :   long tx = typ(x), ty = typ(y);
    2298           0 :   if (tx == ty)
    2299           0 :     return tx == t_INT? subii(x,y): ZC_sub(x,y);
    2300             :   else
    2301           0 :     return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
    2302             : }
    2303             : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
    2304             : static GEN
    2305       27321 : zkmul(GEN x, GEN y)
    2306             : {
    2307       27321 :   long tx = typ(x), ty = typ(y);
    2308       27321 :   if (ty == t_INT)
    2309       20958 :     return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
    2310             :   else
    2311        6363 :     return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
    2312             : }
    2313             : 
    2314             : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
    2315             :  * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
    2316             :  * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
    2317             :  * shallow */
    2318             : GEN
    2319           0 : zkchinese(GEN zkc, GEN x, GEN y)
    2320             : {
    2321           0 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
    2322           0 :   return zk_modHNF(z, UV);
    2323             : }
    2324             : /* special case z = x mod U, = 1 mod V; shallow */
    2325             : GEN
    2326       27321 : zkchinese1(GEN zkc, GEN x)
    2327             : {
    2328       27321 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
    2329       27321 :   return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
    2330             : }
    2331             : static GEN
    2332       25186 : zkVchinese1(GEN zkc, GEN v)
    2333             : {
    2334             :   long i, ly;
    2335       25186 :   GEN y = cgetg_copy(v, &ly);
    2336       25186 :   for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
    2337       25186 :   return y;
    2338             : }
    2339             : 
    2340             : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
    2341             : GEN
    2342       24927 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
    2343             : {
    2344             :   GEN v;
    2345             :   long e;
    2346       24927 :   nf = checknf(nf);
    2347       24927 :   v = idealaddtoone_raw(nf, A, B);
    2348       24927 :   if ((e = gexpo(v)) > 5)
    2349             :   {
    2350        1407 :     GEN b = (typ(v) == t_COL)? v: scalarcol_shallow(v, nf_get_degree(nf));
    2351        1407 :     b= ZC_reducemodlll(b, AB);
    2352        1407 :     if (gexpo(b) < e) v = b;
    2353             :   }
    2354       24927 :   return mkvec2(zk_scalar_or_multable(nf,v), AB);
    2355             : }
    2356             : /* prepare to solve z = x (mod A), z = 1 mod (B)
    2357             :  * and then         z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
    2358             : static GEN
    2359         259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
    2360             : {
    2361         259 :   GEN zkc = zkchineseinit(nf, A, B, AB);
    2362         259 :   GEN mv = gel(zkc,1), mu;
    2363         259 :   if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
    2364          35 :   mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
    2365          35 :   return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
    2366             : }
    2367             : 
    2368             : static GEN
    2369      362306 : apply_U(GEN L, GEN a)
    2370             : {
    2371      362306 :   GEN e, U = gel(L,3), dU = gel(L,4);
    2372      362306 :   if (typ(a) == t_INT)
    2373      126203 :     e = ZC_Z_mul(gel(U,1), subiu(a, 1));
    2374             :   else
    2375             :   { /* t_COL */
    2376      236103 :     GEN t = gel(a,1);
    2377      236103 :     gel(a,1) = subiu(gel(a,1), 1); /* a -= 1 */
    2378      236103 :     e = ZM_ZC_mul(U, a);
    2379      236103 :     gel(a,1) = t; /* restore */
    2380             :   }
    2381      362306 :   return gdiv(e, dU);
    2382             : }
    2383             : 
    2384             : /* true nf; vectors of [[cyc],[g],U.X^-1]. Assume k > 1. */
    2385             : static GEN
    2386       20629 : principal_units(GEN nf, GEN pr, long k, GEN prk)
    2387             : {
    2388             :   GEN list, prb;
    2389       20629 :   ulong mask = quadratic_prec_mask(k);
    2390       20629 :   long a = 1;
    2391             : 
    2392       20629 :   prb = pr_hnf(nf,pr);
    2393       20629 :   list = vectrunc_init(k);
    2394       71155 :   while (mask > 1)
    2395             :   {
    2396       29897 :     GEN pra = prb;
    2397       29897 :     long b = a << 1;
    2398             : 
    2399       29897 :     if (mask & 1) b--;
    2400       29897 :     mask >>= 1;
    2401             :     /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
    2402       29897 :     prb = (b >= k)? prk: idealpows(nf,pr,b);
    2403       29897 :     vectrunc_append(list, zidealij(pra, prb));
    2404       29897 :     a = b;
    2405             :   }
    2406       20629 :   return list;
    2407             : }
    2408             : /* a = 1 mod (pr) return log(a) on local-gens of 1+pr/1+pr^k */
    2409             : static GEN
    2410      228172 : log_prk1(GEN nf, GEN a, long nh, GEN L2, GEN prk)
    2411             : {
    2412      228172 :   GEN y = cgetg(nh+1, t_COL);
    2413      228172 :   long j, iy, c = lg(L2)-1;
    2414      590471 :   for (j = iy = 1; j <= c; j++)
    2415             :   {
    2416      362306 :     GEN L = gel(L2,j), cyc = gel(L,1), gen = gel(L,2), E = apply_U(L,a);
    2417      362306 :     long i, nc = lg(cyc)-1;
    2418      362306 :     int last = (j == c);
    2419     1282064 :     for (i = 1; i <= nc; i++, iy++)
    2420             :     {
    2421      919765 :       GEN t, e = gel(E,i);
    2422      919765 :       if (typ(e) != t_INT) pari_err_COPRIME("zlog_prk1", a, prk);
    2423      919758 :       t = Fp_neg(e, gel(cyc,i));
    2424      919758 :       gel(y,iy) = negi(t);
    2425      919758 :       if (!last && signe(t)) a = nfmulpowmodideal(nf, a, gel(gen,i), t, prk);
    2426             :     }
    2427             :   }
    2428      228165 :   return y;
    2429             : }
    2430             : /* true nf */
    2431             : static GEN
    2432        8057 : principal_units_relations(GEN nf, GEN L2, GEN prk, long nh)
    2433             : {
    2434        8057 :   GEN h = cgetg(nh+1,t_MAT);
    2435        8057 :   long ih, j, c = lg(L2)-1;
    2436       25382 :   for (j = ih = 1; j <= c; j++)
    2437             :   {
    2438       17325 :     GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
    2439       17325 :     long k, lG = lg(G);
    2440       62643 :     for (k = 1; k < lG; k++,ih++)
    2441             :     { /* log(g^f) mod pr^e */
    2442       45318 :       GEN a = nfpowmodideal(nf,gel(G,k),gel(F,k),prk);
    2443       45318 :       gel(h,ih) = ZC_neg(log_prk1(nf, a, nh, L2, prk));
    2444       45318 :       gcoeff(h,ih,ih) = gel(F,k);
    2445             :     }
    2446             :   }
    2447        8057 :   return h;
    2448             : }
    2449             : /* true nf; e > 1; multiplicative group (1 + pr) / (1 + pr^k),
    2450             :  * prk = pr^k or NULL */
    2451             : static GEN
    2452       20629 : idealprincipalunits_i(GEN nf, GEN pr, long k, GEN *pU)
    2453             : {
    2454       20629 :   GEN cyc, gen, L2, prk = idealpows(nf, pr, k);
    2455             : 
    2456       20629 :   L2 = principal_units(nf, pr, k, prk);
    2457       20629 :   if (k == 2)
    2458             :   {
    2459       12572 :     GEN L = gel(L2,1);
    2460       12572 :     cyc = gel(L,1);
    2461       12572 :     gen = gel(L,2);
    2462       12572 :     if (pU) *pU = matid(lg(gen)-1);
    2463             :   }
    2464             :   else
    2465             :   {
    2466        8057 :     long c = lg(L2), j;
    2467        8057 :     GEN EX, h, Ui, vg = cgetg(c, t_VEC);
    2468        8057 :     for (j = 1; j < c; j++) gel(vg, j) = gmael(L2,j,2);
    2469        8057 :     vg = shallowconcat1(vg);
    2470        8057 :     h = principal_units_relations(nf, L2, prk, lg(vg)-1);
    2471        8057 :     h = ZM_hnfall_i(h, NULL, 0);
    2472        8057 :     cyc = ZM_snf_group(h, pU, &Ui);
    2473        8057 :     c = lg(Ui); gen = cgetg(c, t_VEC); EX = gel(cyc,1);
    2474       39844 :     for (j = 1; j < c; j++)
    2475       31787 :       gel(gen,j) = famat_to_nf_modideal_coprime(nf, vg, gel(Ui,j), prk, EX);
    2476             :   }
    2477       20629 :   return mkvec4(cyc, gen, prk, L2);
    2478             : }
    2479             : GEN
    2480         112 : idealprincipalunits(GEN nf, GEN pr, long k)
    2481             : {
    2482             :   pari_sp av;
    2483             :   GEN v;
    2484         112 :   nf = checknf(nf);
    2485         112 :   if (k == 1) { checkprid(pr); retmkvec3(gen_1,cgetg(1,t_VEC),cgetg(1,t_VEC)); }
    2486         105 :   av = avma; v = idealprincipalunits_i(nf, pr, k, NULL);
    2487         105 :   return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), k-1), gel(v,1), gel(v,2)));
    2488             : }
    2489             : 
    2490             : /* true nf; given an ideal pr^k dividing an integral ideal x (in HNF form)
    2491             :  * compute an 'sprk', the structure of G = (Z_K/pr^k)^* [ x = NULL for x=pr^k ]
    2492             :  * Return a vector with at least 4 components [cyc],[gen],[HNF pr^k,pr,k],ff,
    2493             :  * where
    2494             :  * cyc : type of G as abelian group (SNF)
    2495             :  * gen : generators of G, coprime to x
    2496             :  * pr^k: in HNF
    2497             :  * ff  : data for log_g in (Z_K/pr)^*
    2498             :  * Two extra components are present iff k > 1: L2, U
    2499             :  * L2  : list of data structures to compute local DL in (Z_K/pr)^*,
    2500             :  *       and 1 + pr^a/ 1 + pr^b for various a < b <= min(2a, k)
    2501             :  * U   : base change matrices to convert a vector of local DL to DL wrt gen */
    2502             : static GEN
    2503       52934 : sprkinit(GEN nf, GEN pr, long k, GEN x)
    2504             : {
    2505             :   GEN T, p, modpr, cyc, gen, g, g0, ord0, A, prk, U, L2;
    2506       52934 :   long f = pr_get_f(pr);
    2507             : 
    2508       52934 :   if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",k,pr);
    2509       52934 :   modpr = nf_to_Fq_init(nf, &pr,&T,&p);
    2510             :   /* (Z_K / pr)^* */
    2511       52934 :   if (f == 1)
    2512             :   {
    2513       41167 :     g0 = g = pgener_Fp(p);
    2514       41167 :     ord0 = get_arith_ZZM(subiu(p,1));
    2515             :   }
    2516             :   else
    2517             :   {
    2518       11767 :     g0 = g = gener_FpXQ(T,p, &ord0);
    2519       11767 :     g = Fq_to_nf(g, modpr);
    2520       11767 :     if (typ(g) == t_POL) g = poltobasis(nf, g);
    2521             :   }
    2522       52934 :   A = gel(ord0, 1); /* Norm(pr)-1 */
    2523       52934 :   if (k == 1)
    2524             :   {
    2525       32410 :     cyc = mkvec(A);
    2526       32410 :     gen = mkvec(g);
    2527       32410 :     prk = pr_hnf(nf,pr);
    2528       32410 :     L2 = U = NULL;
    2529             :   }
    2530             :   else
    2531             :   { /* local-gens of (1 + pr)/(1 + pr^k) = SNF-gens * U */
    2532             :     GEN AB, B, u, v, w;
    2533             :     long j, l;
    2534       20524 :     w = idealprincipalunits_i(nf, pr, k, &U);
    2535             :     /* incorporate (Z_K/pr)^*, order A coprime to B = expo(1+pr/1+pr^k)*/
    2536       20524 :     cyc = leafcopy(gel(w,1)); B = gel(cyc,1); AB = mulii(A,B);
    2537       20524 :     gen = leafcopy(gel(w,2));
    2538       20524 :     prk = gel(w,3);
    2539       20524 :     g = nfpowmodideal(nf, g, B, prk);
    2540       20524 :     g0 = Fq_pow(g0, modii(B,A), T, p); /* update primitive root */
    2541       20524 :     L2 = mkvec3(A, g, gel(w,4));
    2542       20524 :     gel(cyc,1) = AB;
    2543       20524 :     gel(gen,1) = nfmulmodideal(nf, gel(gen,1), g, prk);
    2544       20524 :     u = mulii(Fp_inv(A,B), A);
    2545       20524 :     v = subui(1, u); l = lg(U);
    2546       20524 :     for (j = 1; j < l; j++) gcoeff(U,1,j) = Fp_mul(u, gcoeff(U,1,j), AB);
    2547       20524 :     U = mkvec2(Rg_col_ei(v, lg(gen)-1, 1), U);
    2548             :   }
    2549             :   /* local-gens of (Z_K/pr^k)^* = SNF-gens * U */
    2550       52934 :   if (x)
    2551             :   {
    2552       24668 :     GEN uv = zkchineseinit(nf, idealmulpowprime(nf,x,pr,utoineg(k)), prk, x);
    2553       24668 :     gen = zkVchinese1(uv, gen);
    2554             :   }
    2555       52934 :   return mkvecn(U? 6: 4, cyc, gen, prk, mkvec3(modpr,g0,ord0), L2, U);
    2556             : }
    2557             : static GEN
    2558      396063 : sprk_get_cyc(GEN s) { return gel(s,1); }
    2559             : static GEN
    2560      123945 : sprk_get_expo(GEN s)
    2561             : {
    2562      123945 :   GEN cyc = sprk_get_cyc(s);
    2563      123945 :   return lg(cyc) == 1? gen_1: gel(cyc, 1);
    2564             : }
    2565             : static GEN
    2566       45038 : sprk_get_gen(GEN s) { return gel(s,2); }
    2567             : static GEN
    2568      306820 : sprk_get_prk(GEN s) { return gel(s,3); }
    2569             : static GEN
    2570      475102 : sprk_get_ff(GEN s) { return gel(s,4); }
    2571             : static GEN
    2572      155984 : sprk_get_pr(GEN s) { GEN ff = gel(s,4); return modpr_get_pr(gel(ff,1)); }
    2573             : /* A = Npr-1, <g> = (Z_K/pr)^*, L2 to 1 + pr / 1 + pr^k */
    2574             : static void
    2575      200823 : sprk_get_L2(GEN s, GEN *A, GEN *g, GEN *L2)
    2576      200823 : { GEN v = gel(s,5); *A = gel(v,1); *g = gel(v,2); *L2 = gel(v,3); }
    2577             : static void
    2578      182854 : sprk_get_U2(GEN s, GEN *U1, GEN *U2)
    2579      182854 : { GEN v = gel(s,6); *U1 = gel(v,1); *U2 = gel(v,2); }
    2580             : static int
    2581      475102 : sprk_is_prime(GEN s) { return lg(s) == 5; }
    2582             : 
    2583             : static GEN
    2584      123903 : famat_zlog_pr(GEN nf, GEN g, GEN e, GEN sprk)
    2585             : {
    2586      123903 :   GEN x = famat_makecoprime(nf, g, e, sprk_get_pr(sprk), sprk_get_prk(sprk),
    2587             :                             sprk_get_expo(sprk));
    2588      123903 :   return log_prk(nf, x, sprk);
    2589             : }
    2590             : /* famat_zlog_pr assuming (g,sprk.pr) = 1 */
    2591             : static GEN
    2592          42 : famat_zlog_pr_coprime(GEN nf, GEN g, GEN e, GEN sprk)
    2593             : {
    2594          42 :   GEN x = famat_to_nf_modideal_coprime(nf, g, e, sprk_get_prk(sprk),
    2595             :                                        sprk_get_expo(sprk));
    2596          42 :   return log_prk(nf, x, sprk);
    2597             : }
    2598             : 
    2599             : /* log_g(a) in (Z_K/pr)^* */
    2600             : static GEN
    2601      475102 : nf_log(GEN nf, GEN a, GEN ff)
    2602             : {
    2603      475102 :   GEN pr = gel(ff,1), g = gel(ff,2), ord = gel(ff,3);
    2604      475102 :   GEN T,p, modpr = nf_to_Fq_init(nf, &pr, &T, &p);
    2605      475102 :   return Fq_log(nf_to_Fq(nf,a,modpr), g, ord, T, p);
    2606             : }
    2607             : /* a in Z_K (t_COL or t_INT), pr prime ideal, sprk = sprkinit(nf,pr,k,x).
    2608             :  * return log(a) on SNF-generators of (Z_K/pr^k)^**/
    2609             : GEN
    2610      479967 : log_prk(GEN nf, GEN a, GEN sprk)
    2611             : {
    2612             :   GEN e, prk, A, g, L2, U1, U2, y;
    2613             : 
    2614      479967 :   if (typ(a) == t_MAT) return famat_zlog_pr(nf, gel(a,1), gel(a,2), sprk);
    2615             : 
    2616      475102 :   e = nf_log(nf, a, sprk_get_ff(sprk));
    2617      475102 :   if (sprk_is_prime(sprk)) return mkcol(e); /* k = 1 */
    2618      182854 :   prk = sprk_get_prk(sprk);
    2619      182854 :   sprk_get_L2(sprk, &A,&g,&L2);
    2620      182854 :   if (signe(e))
    2621             :   {
    2622       52382 :     e = Fp_neg(e, A);
    2623       52382 :     a = nfmulpowmodideal(nf, a, g, e, prk);
    2624             :   }
    2625      182854 :   sprk_get_U2(sprk, &U1,&U2);
    2626      182854 :   y = ZM_ZC_mul(U2, log_prk1(nf, a, lg(U2)-1, L2, prk));
    2627      182847 :   if (signe(e)) y = ZC_sub(y, ZC_Z_mul(U1,e));
    2628      182847 :   return vecmodii(y, sprk_get_cyc(sprk));
    2629             : }
    2630             : GEN
    2631        7896 : log_prk_init(GEN nf, GEN pr, long k)
    2632        7896 : { return sprkinit(checknf(nf),pr,k,NULL);}
    2633             : GEN
    2634         497 : veclog_prk(GEN nf, GEN v, GEN sprk)
    2635             : {
    2636         497 :   long l = lg(v), i;
    2637         497 :   GEN w = cgetg(l, t_MAT);
    2638         497 :   for (i = 1; i < l; i++) gel(w,i) = log_prk(nf, gel(v,i), sprk);
    2639         497 :   return w;
    2640             : }
    2641             : 
    2642             : static GEN
    2643      118550 : famat_zlog(GEN nf, GEN fa, GEN sgn, zlog_S *S)
    2644             : {
    2645      118550 :   long i, l0, l = lg(S->U);
    2646      118550 :   GEN g = gel(fa,1), e = gel(fa,2), y = cgetg(l, t_COL);
    2647      118550 :   l0 = lg(S->sprk); /* = l (trivial arch. part), or l-1 */
    2648      118550 :   for (i=1; i < l0; i++) gel(y,i) = famat_zlog_pr(nf, g, e, gel(S->sprk,i));
    2649      118550 :   if (l0 != l)
    2650             :   {
    2651       93566 :     if (!sgn) sgn = nfsign_arch(nf, fa, S->archp);
    2652       93566 :     gel(y,l0) = Flc_to_ZC(sgn);
    2653             :   }
    2654      118550 :   return y;
    2655             : }
    2656             : 
    2657             : /* assume that cyclic factors are normalized, in particular != [1] */
    2658             : static GEN
    2659       44660 : split_U(GEN U, GEN Sprk)
    2660             : {
    2661       44660 :   long t = 0, k, n, l = lg(Sprk);
    2662       44660 :   GEN vU = cgetg(l+1, t_VEC);
    2663       88893 :   for (k = 1; k < l; k++)
    2664             :   {
    2665       44233 :     n = lg(sprk_get_cyc(gel(Sprk,k))) - 1; /* > 0 */
    2666       44233 :     gel(vU,k) = vecslice(U, t+1, t+n);
    2667       44233 :     t += n;
    2668             :   }
    2669             :   /* t+1 .. lg(U)-1 */
    2670       44660 :   n = lg(U) - t - 1; /* can be 0 */
    2671       44660 :   if (!n) setlg(vU,l); else gel(vU,l) = vecslice(U, t+1, t+n);
    2672       44660 :   return vU;
    2673             : }
    2674             : 
    2675             : void
    2676      432906 : init_zlog(zlog_S *S, GEN bid)
    2677             : {
    2678      432906 :   GEN fa2 = bid_get_fact2(bid);
    2679      432906 :   S->U = bid_get_U(bid);
    2680      432906 :   S->hU = lg(bid_get_cyc(bid))-1;
    2681      432906 :   S->archp = bid_get_archp(bid);
    2682      432906 :   S->sprk = bid_get_sprk(bid);
    2683      432906 :   S->bid = bid;
    2684      432906 :   S->P = gel(fa2,1);
    2685      432906 :   S->k = gel(fa2,2);
    2686      432906 :   S->no2 = lg(S->P) == lg(gel(bid_get_fact(bid),1));
    2687      432906 : }
    2688             : 
    2689             : /* a a t_FRAC/t_INT, reduce mod bid */
    2690             : static GEN
    2691          14 : Q_mod_bid(GEN bid, GEN a)
    2692             : {
    2693          14 :   GEN xZ = gcoeff(bid_get_ideal(bid),1,1);
    2694          14 :   GEN b = Rg_to_Fp(a, xZ);
    2695          14 :   if (gsigne(a) < 0) b = subii(b, xZ);
    2696          14 :   return signe(b)? b: xZ;
    2697             : }
    2698             : /* Return decomposition of a on the CRT generators blocks attached to the
    2699             :  * S->sprk and sarch; sgn = sign(a, S->arch), NULL if unknown */
    2700             : static GEN
    2701      306191 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S)
    2702             : {
    2703             :   long k, l;
    2704             :   GEN y;
    2705      306191 :   a = nf_to_scalar_or_basis(nf, a);
    2706      306191 :   switch(typ(a))
    2707             :   {
    2708       85911 :     case t_INT: break;
    2709          14 :     case t_FRAC: a = Q_mod_bid(S->bid, a); break;
    2710             :     default: /* case t_COL: */
    2711             :     {
    2712             :       GEN den;
    2713      220266 :       check_nfelt(a, &den);
    2714      220266 :       if (den)
    2715             :       {
    2716       68597 :         a = Q_muli_to_int(a, den);
    2717       68597 :         a = mkmat2(mkcol2(a, den), mkcol2(gen_1, gen_m1));
    2718       68597 :         return famat_zlog(nf, a, sgn, S);
    2719             :       }
    2720             :     }
    2721             :   }
    2722      237594 :   if (sgn)
    2723       60669 :     sgn = (lg(sgn) == 1)? NULL: leafcopy(sgn);
    2724             :   else
    2725      176925 :     sgn = (lg(S->archp) == 1)? NULL: nfsign_arch(nf, a, S->archp);
    2726      237594 :   l = lg(S->sprk);
    2727      237594 :   y = cgetg(sgn? l+1: l, t_COL);
    2728      536930 :   for (k = 1; k < l; k++)
    2729             :   {
    2730      299343 :     GEN sprk = gel(S->sprk,k);
    2731      299343 :     gel(y,k) = log_prk(nf, a, sprk);
    2732             :   }
    2733      237587 :   if (sgn) gel(y,l) = Flc_to_ZC(sgn);
    2734      237587 :   return y;
    2735             : }
    2736             : 
    2737             : /* true nf */
    2738             : GEN
    2739       14343 : pr_basis_perm(GEN nf, GEN pr)
    2740             : {
    2741       14343 :   long f = pr_get_f(pr);
    2742             :   GEN perm;
    2743       14343 :   if (f == nf_get_degree(nf)) return identity_perm(f);
    2744       12614 :   perm = cgetg(f+1, t_VECSMALL);
    2745       12614 :   perm[1] = 1;
    2746       12614 :   if (f > 1)
    2747             :   {
    2748         539 :     GEN H = pr_hnf(nf,pr);
    2749             :     long i, k;
    2750        1743 :     for (i = k = 2; k <= f; i++)
    2751        1204 :       if (!equali1(gcoeff(H,i,i))) perm[k++] = i;
    2752             :   }
    2753       12614 :   return perm;
    2754             : }
    2755             : 
    2756             : /* \sum U[i]*y[i], U[i] ZM, y[i] ZC. We allow lg(y) > lg(U). */
    2757             : static GEN
    2758      356137 : ZMV_ZCV_mul(GEN U, GEN y)
    2759             : {
    2760      356137 :   long i, l = lg(U);
    2761      356137 :   GEN z = NULL;
    2762      356137 :   if (l == 1) return cgetg(1,t_COL);
    2763      992167 :   for (i = 1; i < l; i++)
    2764             :   {
    2765      636030 :     GEN u = ZM_ZC_mul(gel(U,i), gel(y,i));
    2766      636030 :     z = z? ZC_add(z, u): u;
    2767             :   }
    2768      356137 :   return z;
    2769             : }
    2770             : /* A * (U[1], ..., U[d] */
    2771             : static GEN
    2772         518 : ZM_ZMV_mul(GEN A, GEN U)
    2773             : {
    2774             :   long i, l;
    2775         518 :   GEN V = cgetg_copy(U,&l);
    2776         518 :   for (i = 1; i < l; i++) gel(V,i) = ZM_mul(A,gel(U,i));
    2777         518 :   return V;
    2778             : }
    2779             : 
    2780             : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
    2781             :  * defined implicitly via CRT. 'ind' is the index of pr in modulus
    2782             :  * factorization */
    2783             : GEN
    2784       93471 : log_gen_pr(zlog_S *S, long ind, GEN nf, long e)
    2785             : {
    2786       93471 :   GEN A, sprk = gel(S->sprk,ind), Uind = gel(S->U, ind);
    2787             : 
    2788       93471 :   if (e == 1) retmkmat( gel(Uind,1) );
    2789             :   else
    2790             :   {
    2791       32060 :     GEN G, pr = sprk_get_pr(sprk);
    2792             :     long i, l;
    2793       32060 :     if (e == 2)
    2794             :     {
    2795       17969 :       GEN A, g, L, L2; sprk_get_L2(sprk,&A,&g,&L2); L = gel(L2,1);
    2796       17969 :       G = gel(L,2); l = lg(G);
    2797             :     }
    2798             :     else
    2799             :     {
    2800       14091 :       GEN perm = pr_basis_perm(nf,pr), PI = nfpow_u(nf, pr_get_gen(pr), e-1);
    2801       14091 :       l = lg(perm);
    2802       14091 :       G = cgetg(l, t_VEC);
    2803       14091 :       if (typ(PI) == t_INT)
    2804             :       { /* zk_ei_mul doesn't allow t_INT */
    2805        1722 :         long N = nf_get_degree(nf);
    2806        1722 :         gel(G,1) = addiu(PI,1);
    2807        2898 :         for (i = 2; i < l; i++)
    2808             :         {
    2809        1176 :           GEN z = col_ei(N, 1);
    2810        1176 :           gel(G,i) = z; gel(z, perm[i]) = PI;
    2811             :         }
    2812             :       }
    2813             :       else
    2814             :       {
    2815       12369 :         gel(G,1) = nfadd(nf, gen_1, PI);
    2816       12719 :         for (i = 2; i < l; i++)
    2817         350 :           gel(G,i) = nfadd(nf, gen_1, zk_ei_mul(nf, PI, perm[i]));
    2818             :       }
    2819             :     }
    2820       32060 :     A = cgetg(l, t_MAT);
    2821       68327 :     for (i = 1; i < l; i++)
    2822       36267 :       gel(A,i) = ZM_ZC_mul(Uind, log_prk(nf, gel(G,i), sprk));
    2823       32060 :     return A;
    2824             :   }
    2825             : }
    2826             : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
    2827             :  * v = vector of r1 real places */
    2828             : GEN
    2829       15043 : log_gen_arch(zlog_S *S, long index)
    2830             : {
    2831       15043 :   GEN U = gel(S->U, lg(S->U)-1);
    2832       15043 :   return gel(U, index);
    2833             : }
    2834             : 
    2835             : /* compute bid.clgp: [h,cyc] or [h,cyc,gen] */
    2836             : static GEN
    2837       45731 : bid_grp(GEN nf, GEN U, GEN cyc, GEN g, GEN F, GEN sarch)
    2838             : {
    2839       45731 :   GEN G, h = ZV_prod(cyc);
    2840             :   long c;
    2841       45731 :   if (!U) return mkvec2(h,cyc);
    2842       45472 :   c = lg(U);
    2843       45472 :   G = cgetg(c,t_VEC);
    2844       45472 :   if (c > 1)
    2845             :   {
    2846       38479 :     GEN U0, Uoo, EX = gel(cyc,1); /* exponent of bid */
    2847       38479 :     long i, hU = nbrows(U), nba = lg(sarch_get_cyc(sarch))-1; /* #f_oo */
    2848       38479 :     if (!nba) { U0 = U; Uoo = NULL; }
    2849       20790 :     else if (nba == hU) { U0 = NULL; Uoo = U; }
    2850             :     else
    2851             :     {
    2852       16989 :       U0 = rowslice(U, 1, hU-nba);
    2853       16989 :       Uoo = rowslice(U, hU-nba+1, hU);
    2854             :     }
    2855      111076 :     for (i = 1; i < c; i++)
    2856             :     {
    2857       72597 :       GEN t = gen_1;
    2858       72597 :       if (U0) t = famat_to_nf_modideal_coprime(nf, g, gel(U0,i), F, EX);
    2859       72597 :       if (Uoo) t = set_sign_mod_divisor(nf, ZV_to_Flv(gel(Uoo,i),2), t, sarch);
    2860       72597 :       gel(G,i) = t;
    2861             :     }
    2862             :   }
    2863       45472 :   return mkvec3(h, cyc, G);
    2864             : }
    2865             : 
    2866             : /* remove prime ideals of norm 2 with exponent 1 from factorization */
    2867             : static GEN
    2868       45430 : famat_strip2(GEN fa)
    2869             : {
    2870       45430 :   GEN P = gel(fa,1), E = gel(fa,2), Q, F;
    2871       45430 :   long l = lg(P), i, j;
    2872       45430 :   Q = cgetg(l, t_COL);
    2873       45430 :   F = cgetg(l, t_COL);
    2874       98497 :   for (i = j = 1; i < l; i++)
    2875             :   {
    2876       53067 :     GEN pr = gel(P,i), e = gel(E,i);
    2877       53067 :     if (!absequaliu(pr_get_p(pr), 2) || itou(e) != 1 || pr_get_f(pr) != 1)
    2878             :     {
    2879       45038 :       gel(Q,j) = pr;
    2880       45038 :       gel(F,j) = e; j++;
    2881             :     }
    2882             :   }
    2883       45430 :   setlg(Q,j);
    2884       45430 :   setlg(F,j); return mkmat2(Q,F);
    2885             : }
    2886             : 
    2887             : /* Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
    2888             :    flag may include nf_GEN | nf_INIT */
    2889             : static GEN
    2890       45451 : Idealstar_i(GEN nf, GEN ideal, long flag)
    2891             : {
    2892             :   long i, k, nbp, R1;
    2893       45451 :   GEN y, cyc, U, u1 = NULL, fa, fa2, sprk, x, arch, archp, E, P, sarch, gen;
    2894             : 
    2895       45451 :   nf = checknf(nf);
    2896       45451 :   R1 = nf_get_r1(nf);
    2897       45451 :   if (typ(ideal) == t_VEC && lg(ideal) == 3)
    2898             :   {
    2899       22456 :     arch = gel(ideal,2);
    2900       22456 :     ideal= gel(ideal,1);
    2901       22456 :     switch(typ(arch))
    2902             :     {
    2903             :       case t_VEC:
    2904       22008 :         if (lg(arch) != R1+1)
    2905           0 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    2906       22008 :         archp = vec01_to_indices(arch);
    2907       22008 :         break;
    2908             :       case t_VECSMALL:
    2909         448 :         archp = arch;
    2910         448 :         k = lg(archp)-1;
    2911         448 :         if (k && archp[k] > R1)
    2912           7 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    2913         441 :         arch = indices_to_vec01(archp, R1);
    2914         441 :         break;
    2915             :       default:
    2916           0 :         pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    2917             :         return NULL;/*LCOV_EXCL_LINE*/
    2918             :     }
    2919       22449 :   }
    2920             :   else
    2921             :   {
    2922       22995 :     arch = zerovec(R1);
    2923       22995 :     archp = cgetg(1, t_VECSMALL);
    2924             :   }
    2925       45444 :   if (is_nf_factor(ideal))
    2926             :   {
    2927        1148 :     fa = ideal;
    2928        1148 :     x = idealfactorback(nf, gel(fa,1), gel(fa,2), 0);
    2929             :   }
    2930             :   else
    2931             :   {
    2932       44296 :     fa = idealfactor(nf, ideal);
    2933       44289 :     x = ideal;
    2934             :   }
    2935       45437 :   if (typ(x) != t_MAT) x = idealhnf_shallow(nf, x);
    2936       45437 :   if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
    2937       45437 :   if (typ(gcoeff(x,1,1)) != t_INT)
    2938           7 :     pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
    2939       45430 :   sarch = nfarchstar(nf, x, archp);
    2940       45430 :   fa2 = famat_strip2(fa);
    2941       45430 :   P = gel(fa2,1);
    2942       45430 :   E = gel(fa2,2);
    2943       45430 :   nbp = lg(P)-1;
    2944       45430 :   sprk = cgetg(nbp+1,t_VEC);
    2945       45430 :   if (nbp)
    2946             :   {
    2947       34608 :     GEN t = (lg(gel(fa,1))==2)? NULL: x; /* beware fa != fa2 */
    2948       34608 :     cyc = cgetg(nbp+2,t_VEC);
    2949       34608 :     gen = cgetg(nbp+1,t_VEC);
    2950       79646 :     for (i = 1; i <= nbp; i++)
    2951             :     {
    2952       45038 :       GEN L = sprkinit(nf, gel(P,i), itou(gel(E,i)), t);
    2953       45038 :       gel(sprk,i) = L;
    2954       45038 :       gel(cyc,i) = sprk_get_cyc(L);
    2955             :       /* true gens are congruent to those mod x AND positive at archp */
    2956       45038 :       gel(gen,i) = sprk_get_gen(L);
    2957             :     }
    2958       34608 :     gel(cyc,i) = sarch_get_cyc(sarch);
    2959       34608 :     cyc = shallowconcat1(cyc);
    2960       34608 :     gen = shallowconcat1(gen);
    2961       34608 :     cyc = ZV_snf_group(cyc, &U, (flag & nf_GEN)? &u1: NULL);
    2962             :   }
    2963             :   else
    2964             :   {
    2965       10822 :     cyc = sarch_get_cyc(sarch);
    2966       10822 :     gen = cgetg(1,t_VEC);
    2967       10822 :     U = matid(lg(cyc)-1);
    2968       10822 :     if (flag & nf_GEN) u1 = U;
    2969             :   }
    2970       45430 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    2971       45430 :   if (!(flag & nf_INIT)) return y;
    2972       44618 :   U = split_U(U, sprk);
    2973       44618 :   return mkvec5(mkvec2(x, arch), y, mkvec2(fa,fa2), mkvec2(sprk, sarch), U);
    2974             : }
    2975             : GEN
    2976       45192 : Idealstar(GEN nf, GEN ideal, long flag)
    2977             : {
    2978       45192 :   pari_sp av = avma;
    2979       45192 :   if (!nf) nf = nfinit(pol_x(0), DEFAULTPREC);
    2980       45192 :   return gerepilecopy(av, Idealstar_i(nf, ideal, flag));
    2981             : }
    2982             : GEN
    2983         259 : Idealstarprk(GEN nf, GEN pr, long k, long flag)
    2984             : {
    2985         259 :   pari_sp av = avma;
    2986         259 :   GEN z = Idealstar_i(nf, mkmat2(mkcol(pr),mkcols(k)), flag);
    2987         259 :   return gerepilecopy(av, z);
    2988             : }
    2989             : 
    2990             : /* FIXME: obsolete */
    2991             : GEN
    2992           0 : zidealstarinitgen(GEN nf, GEN ideal)
    2993           0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
    2994             : GEN
    2995           0 : zidealstarinit(GEN nf, GEN ideal)
    2996           0 : { return Idealstar(nf,ideal, nf_INIT); }
    2997             : GEN
    2998           0 : zidealstar(GEN nf, GEN ideal)
    2999           0 : { return Idealstar(nf,ideal, nf_GEN); }
    3000             : 
    3001             : GEN
    3002          70 : idealstar0(GEN nf, GEN ideal,long flag)
    3003             : {
    3004          70 :   switch(flag)
    3005             :   {
    3006           0 :     case 0: return Idealstar(nf,ideal, nf_GEN);
    3007          56 :     case 1: return Idealstar(nf,ideal, nf_INIT);
    3008          14 :     case 2: return Idealstar(nf,ideal, nf_INIT|nf_GEN);
    3009           0 :     default: pari_err_FLAG("idealstar");
    3010             :   }
    3011             :   return NULL; /* LCOV_EXCL_LINE */
    3012             : }
    3013             : 
    3014             : void
    3015      220266 : check_nfelt(GEN x, GEN *den)
    3016             : {
    3017      220266 :   long l = lg(x), i;
    3018      220266 :   GEN t, d = NULL;
    3019      220266 :   if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
    3020      814262 :   for (i=1; i<l; i++)
    3021             :   {
    3022      593996 :     t = gel(x,i);
    3023      593996 :     switch (typ(t))
    3024             :     {
    3025      448822 :       case t_INT: break;
    3026             :       case t_FRAC:
    3027      145174 :         if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
    3028      145174 :         break;
    3029           0 :       default: pari_err_TYPE("check_nfelt", x);
    3030             :     }
    3031             :   }
    3032      220266 :   *den = d;
    3033      220266 : }
    3034             : 
    3035             : GEN
    3036     1607665 : vecmodii(GEN x, GEN y)
    3037     1607665 : { pari_APPLY_same(modii(gel(x,i), gel(y,i))) }
    3038             : 
    3039             : GEN
    3040       27230 : vecmoduu(GEN x, GEN y)
    3041       27230 : { pari_APPLY_ulong(uel(x,i) % uel(y,i)) }
    3042             : 
    3043             : /* assume a true bnf and bid */
    3044             : GEN
    3045       37429 : ideallog_units(GEN bnf, GEN bid)
    3046             : {
    3047       37429 :   GEN nf = bnf_get_nf(bnf), D, y, C, cyc;
    3048       37429 :   long j, lU = lg(bnf_get_logfu(bnf)); /* r1+r2 */
    3049             :   zlog_S S;
    3050       37429 :   init_zlog(&S, bid);
    3051       37429 :   if (!S.hU) return zeromat(0,lU);
    3052       37429 :   cyc = bid_get_cyc(bid);
    3053       37429 :   D = nfsign_fu(bnf, bid_get_archp(bid));
    3054       37429 :   y = cgetg(lU, t_MAT);
    3055       37429 :   if ((C = bnf_build_cheapfu(bnf)))
    3056       37429 :   { for (j = 1; j < lU; j++) gel(y,j) = zlog(nf, gel(C,j), gel(D,j), &S); }
    3057             :   else
    3058             :   {
    3059           0 :     long i, l = lg(S.U), l0 = lg(S.sprk);
    3060             :     GEN X, U;
    3061           0 :     if (!(C = bnf_compactfu_mat(bnf))) bnf_build_units(bnf); /* error */
    3062           0 :     X = gel(C,1); U = gel(C,2);
    3063           0 :     for (j = 1; j < lU; j++) gel(y,j) = cgetg(l, t_COL);
    3064           0 :     for (i = 1; i < l0; i++)
    3065             :     {
    3066           0 :       GEN sprk = gel(S.sprk, i);
    3067           0 :       GEN Xi = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3068           0 :       for (j = 1; j < lU; j++)
    3069           0 :         gcoeff(y,i,j) = famat_zlog_pr_coprime(nf, Xi, gel(U,j), sprk);
    3070             :     }
    3071           0 :     if (l0 != l)
    3072           0 :       for (j = 1; j < lU; j++) gcoeff(y,l0,j) = Flc_to_ZC(gel(D,j));
    3073             :   }
    3074       37429 :   y = vec_prepend(y, zlog(nf, bnf_get_tuU(bnf), nfsign_tu(bnf, S.archp), &S));
    3075       98098 :   for (j = 1; j <= lU; j++)
    3076       60669 :     gel(y,j) = vecmodii(ZMV_ZCV_mul(S.U, gel(y,j)), cyc);
    3077       37429 :   return y;
    3078             : }
    3079             : GEN
    3080         518 : log_prk_units(GEN nf, GEN D, GEN sprk)
    3081             : {
    3082         518 :   GEN L, Ltu = log_prk(nf, gel(D,1), sprk);
    3083         518 :   D = gel(D,2);
    3084         518 :   if (lg(D) != 3 || typ(gel(D,2)) != t_MAT) L = veclog_prk(nf, D, sprk);
    3085             :   else
    3086             :   {
    3087          21 :     GEN X = gel(D,1), U = gel(D,2);
    3088          21 :     long j, lU = lg(U);
    3089          21 :     X = sunits_makecoprime(X, sprk_get_pr(sprk), sprk_get_prk(sprk));
    3090          21 :     L = cgetg(lU, t_MAT);
    3091          63 :     for (j = 1; j < lU; j++)
    3092          42 :       gel(L,j) = famat_zlog_pr_coprime(nf, X, gel(U,j), sprk);
    3093             :   }
    3094         518 :   return vec_prepend(L, Ltu);
    3095             : }
    3096             : 
    3097             : static GEN
    3098      297218 : ideallog_i(GEN nf, GEN x, zlog_S *S)
    3099             : {
    3100      297218 :   pari_sp av = avma;
    3101             :   GEN y;
    3102      297218 :   if (!S->hU) return cgetg(1, t_COL);
    3103      295475 :   if (typ(x) == t_MAT)
    3104       49953 :     y = famat_zlog(nf, x, NULL, S);
    3105             :   else
    3106      245522 :     y = zlog(nf, x, NULL, S);
    3107      295468 :   y = ZMV_ZCV_mul(S->U, y);
    3108      295468 :   return gerepileupto(av, vecmodii(y, bid_get_cyc(S->bid)));
    3109             : }
    3110             : GEN
    3111      303896 : ideallog(GEN nf, GEN x, GEN bid)
    3112             : {
    3113             :   zlog_S S;
    3114      303896 :   if (!nf) return Zideallog(bid, x);
    3115      297225 :   checkbid(bid); init_zlog(&S, bid);
    3116      297218 :   return ideallog_i(checknf(nf), x, &S);
    3117             : }
    3118             : 
    3119             : /*************************************************************************/
    3120             : /**                                                                     **/
    3121             : /**               JOIN BID STRUCTURES, IDEAL LISTS                      **/
    3122             : /**                                                                     **/
    3123             : /*************************************************************************/
    3124             : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
    3125             :  * Output: bid for m1 m2 */
    3126             : static GEN
    3127         462 : join_bid(GEN nf, GEN bid1, GEN bid2)
    3128             : {
    3129         462 :   pari_sp av = avma;
    3130             :   long nbgen, l1,l2;
    3131             :   GEN I1,I2, G1,G2, sprk1,sprk2, cyc1,cyc2, sarch;
    3132         462 :   GEN sprk, fa,fa2, U, cyc, y, u1 = NULL, x, gen;
    3133             : 
    3134         462 :   I1 = bid_get_ideal(bid1);
    3135         462 :   I2 = bid_get_ideal(bid2);
    3136         462 :   if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
    3137         259 :   G1 = bid_get_grp(bid1);
    3138         259 :   G2 = bid_get_grp(bid2);
    3139         259 :   x = idealmul(nf, I1,I2);
    3140         259 :   fa = famat_mul_shallow(bid_get_fact(bid1), bid_get_fact(bid2));
    3141         259 :   fa2= famat_mul_shallow(bid_get_fact2(bid1), bid_get_fact2(bid2));
    3142         259 :   sprk1 = bid_get_sprk(bid1);
    3143         259 :   sprk2 = bid_get_sprk(bid2);
    3144         259 :   sprk = shallowconcat(sprk1, sprk2);
    3145             : 
    3146         259 :   cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
    3147         259 :   cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
    3148         259 :   gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
    3149         259 :   nbgen = l1+l2-2;
    3150         259 :   cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
    3151         259 :   if (nbgen)
    3152             :   {
    3153         259 :     GEN U1 = bid_get_U(bid1), U2 = bid_get_U(bid2);
    3154         259 :     U1 = l1==1? const_vec(lg(sprk1), cgetg(1,t_MAT))
    3155         259 :               : ZM_ZMV_mul(vecslice(U, 1, l1-1),   U1);
    3156         259 :     U2 = l2==1? const_vec(lg(sprk2), cgetg(1,t_MAT))
    3157         259 :               : ZM_ZMV_mul(vecslice(U, l1, nbgen), U2);
    3158         259 :     U = shallowconcat(U1, U2);
    3159             :   }
    3160             :   else
    3161           0 :     U = const_vec(lg(sprk), cgetg(1,t_MAT));
    3162             : 
    3163         259 :   if (gen)
    3164             :   {
    3165         259 :     GEN uv = zkchinese1init2(nf, I2, I1, x);
    3166         518 :     gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
    3167         259 :                         zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
    3168             :   }
    3169         259 :   sarch = bid_get_sarch(bid1); /* trivial */
    3170         259 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3171         259 :   x = mkvec2(x, bid_get_arch(bid1));
    3172         259 :   y = mkvec5(x, y, mkvec2(fa, fa2), mkvec2(sprk, sarch), U);
    3173         259 :   return gerepilecopy(av,y);
    3174             : }
    3175             : 
    3176             : typedef struct _ideal_data {
    3177             :   GEN nf, emb, L, pr, prL, sgnU, archp;
    3178             : } ideal_data;
    3179             : 
    3180             : /* z <- ( z | f(v[i])_{i=1..#v} ) */
    3181             : static void
    3182      119245 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
    3183             : {
    3184      119245 :   long i, nz, lv = lg(v);
    3185             :   GEN z, Z;
    3186      119245 :   if (lv == 1) return;
    3187       52535 :   z = *pz; nz = lg(z)-1;
    3188       52535 :   *pz = Z = cgetg(lv + nz, typ(z));
    3189       52535 :   for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
    3190       52535 :   Z += nz;
    3191       52535 :   for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
    3192             : }
    3193             : static GEN
    3194         462 : join_idealinit(ideal_data *D, GEN x)
    3195         462 : { return join_bid(D->nf, x, D->prL); }
    3196             : static GEN
    3197       60039 : join_ideal(ideal_data *D, GEN x)
    3198       60039 : { return idealmulpowprime(D->nf, x, D->pr, D->L); }
    3199             : static GEN
    3200         441 : join_unit(ideal_data *D, GEN x)
    3201             : {
    3202         441 :   GEN bid = join_idealinit(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3203         441 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3204         441 :   return mkvec2(bid, v);
    3205             : }
    3206             : 
    3207             : GEN
    3208          42 : log_prk_units_init(GEN bnf)
    3209             : {
    3210          42 :   GEN U = bnf_has_fu(bnf);
    3211          42 :   if (U) U = matalgtobasis(bnf_get_nf(bnf), U);
    3212          14 :   else if (!(U = bnf_compactfu_mat(bnf))) (void)bnf_build_units(bnf);
    3213          42 :   return mkvec2(bnf_get_tuU(bnf), U);
    3214             : }
    3215             : /*  flag & nf_GEN : generators, otherwise no
    3216             :  *  flag &2 : units, otherwise no
    3217             :  *  flag &4 : ideals in HNF, otherwise bid
    3218             :  *  flag &8 : omit ideals which cannot be conductors (pr^1 with Npr=2) */
    3219             : static GEN
    3220        6034 : Ideallist(GEN bnf, ulong bound, long flag)
    3221             : {
    3222        6034 :   const long cond = flag & 8;
    3223        6034 :   const long do_units = flag & 2, big_id = !(flag & 4);
    3224        6034 :   const long istar_flag = (flag & nf_GEN) | nf_INIT;
    3225        6034 :   pari_sp av, av0 = avma;
    3226             :   long i, j;
    3227        6034 :   GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
    3228             :   forprime_t S;
    3229             :   ideal_data ID;
    3230        6034 :   GEN (*join_z)(ideal_data*, GEN) =
    3231             :     do_units? &join_unit
    3232        6034 :               : (big_id? &join_idealinit: &join_ideal);
    3233             : 
    3234        6034 :   nf = checknf(bnf);
    3235        6034 :   if ((long)bound <= 0) return empty;
    3236        6034 :   id = matid(nf_get_degree(nf));
    3237        6034 :   if (big_id) id = Idealstar(nf,id, istar_flag);
    3238             : 
    3239             :   /* z[i] will contain all "objects" of norm i. Depending on flag, this means
    3240             :    * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
    3241             :    * in vectors, computed one primary component at a time; join_z
    3242             :    * reconstructs the global object */
    3243        6034 :   BOUND = utoipos(bound);
    3244        6034 :   z = const_vec(bound, empty);
    3245        6034 :   U = do_units? log_prk_units_init(bnf): NULL;
    3246        6034 :   gel(z,1) = mkvec(U? mkvec2(id, cgetg(1,t_VEC)): id);
    3247        6034 :   ID.nf = nf;
    3248             : 
    3249        6034 :   p = cgetipos(3);
    3250        6034 :   u_forprime_init(&S, 2, bound);
    3251        6034 :   av = avma;
    3252       33474 :   while ((p[2] = u_forprime_next(&S)))
    3253             :   {
    3254       21406 :     if (DEBUGLEVEL>1) err_printf("%ld ",p[2]);
    3255       21406 :     fa = idealprimedec_limit_norm(nf, p, BOUND);
    3256       43638 :     for (j = 1; j < lg(fa); j++)
    3257             :     {
    3258       22232 :       GEN pr = gel(fa,j), z2 = leafcopy(z);
    3259       22232 :       ulong Q, q = upr_norm(pr);
    3260             :       long l;
    3261       22232 :       ID.pr = ID.prL = pr;
    3262       22232 :       if (cond && q == 2) { l = 2; Q = 4; } else { l = 1; Q = q; }
    3263       53781 :       for (; Q <= bound; l++, Q *= q) /* add pr^l */
    3264             :       {
    3265             :         ulong iQ;
    3266       31549 :         ID.L = utoipos(l);
    3267       31549 :         if (big_id) {
    3268         203 :           ID.prL = Idealstarprk(nf, pr, l, istar_flag);
    3269         203 :           if (U)
    3270         182 :             ID.emb = Q == 2? cgetg(1,t_VEC)
    3271         182 :                            : log_prk_units(nf, U, gel(bid_get_sprk(ID.prL),1));
    3272             :         }
    3273      150794 :         for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
    3274      119245 :           concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
    3275             :       }
    3276             :     }
    3277       21406 :     if (gc_needed(av,1))
    3278             :     {
    3279           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
    3280           0 :       z = gerepilecopy(av, z);
    3281             :     }
    3282             :   }
    3283        6034 :   return gerepilecopy(av0, z);
    3284             : }
    3285             : GEN
    3286         350 : ideallist0(GEN bnf,long bound, long flag) {
    3287         350 :   if (flag<0 || flag>15) pari_err_FLAG("ideallist");
    3288         350 :   return Ideallist(bnf,bound,flag);
    3289             : }
    3290             : GEN
    3291        5684 : ideallist(GEN bnf,long bound) { return Ideallist(bnf,bound,4); }
    3292             : 
    3293             : /* bid = for module m (without arch. part), arch = archimedean part.
    3294             :  * Output: bid for [m,arch] */
    3295             : static GEN
    3296          42 : join_bid_arch(GEN nf, GEN bid, GEN archp)
    3297             : {
    3298          42 :   pari_sp av = avma;
    3299             :   GEN G, U;
    3300          42 :   GEN sprk, cyc, y, u1 = NULL, x, sarch, gen;
    3301             : 
    3302          42 :   checkbid(bid);
    3303          42 :   G = bid_get_grp(bid);
    3304          42 :   x = bid_get_ideal(bid);
    3305          42 :   sarch = nfarchstar(nf, bid_get_ideal(bid), archp);
    3306          42 :   sprk = bid_get_sprk(bid);
    3307             : 
    3308          42 :   gen = (lg(G)>3)? gel(G,3): NULL;
    3309          42 :   cyc = diagonal_shallow(shallowconcat(gel(G,2), sarch_get_cyc(sarch)));
    3310          42 :   cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
    3311          42 :   y = bid_grp(nf, u1, cyc, gen, x, sarch);
    3312          42 :   U = split_U(U, sprk);
    3313          42 :   y = mkvec5(mkvec2(x, archp), y, gel(bid,3), mkvec2(sprk, sarch), U);
    3314          42 :   return gerepilecopy(av,y);
    3315             : }
    3316             : static GEN
    3317          42 : join_arch(ideal_data *D, GEN x) {
    3318          42 :   return join_bid_arch(D->nf, x, D->archp);
    3319             : }
    3320             : static GEN
    3321          14 : join_archunit(ideal_data *D, GEN x) {
    3322          14 :   GEN bid = join_arch(D, gel(x,1)), u = gel(x,2), v = mkvec(D->emb);
    3323          14 :   if (lg(u) != 1) v = shallowconcat(u, v);
    3324          14 :   return mkvec2(bid, v);
    3325             : }
    3326             : 
    3327             : /* L from ideallist, add archimedean part */
    3328             : GEN
    3329          14 : ideallistarch(GEN bnf, GEN L, GEN arch)
    3330             : {
    3331             :   pari_sp av;
    3332          14 :   long i, j, l = lg(L), lz;
    3333             :   GEN v, z, V;
    3334             :   ideal_data ID;
    3335             :   GEN (*join_z)(ideal_data*, GEN);
    3336             : 
    3337          14 :   if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
    3338          14 :   if (l == 1) return cgetg(1,t_VEC);
    3339          14 :   z = gel(L,1);
    3340          14 :   if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3341          14 :   z = gel(z,1); /* either a bid or [bid,U] */
    3342          14 :   ID.nf = checknf(bnf);
    3343          14 :   ID.archp = vec01_to_indices(arch);
    3344          14 :   if (lg(z) == 3) { /* the latter: do units */
    3345           7 :     if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    3346           7 :     ID.emb = zm_to_ZM( rowpermute(nfsign_units(bnf,NULL,1), ID.archp) );
    3347           7 :     join_z = &join_archunit;
    3348             :   } else
    3349           7 :     join_z = &join_arch;
    3350          14 :   av = avma; V = cgetg(l, t_VEC);
    3351          63 :   for (i = 1; i < l; i++)
    3352             :   {
    3353          49 :     z = gel(L,i); lz = lg(z);
    3354          49 :     gel(V,i) = v = cgetg(lz,t_VEC);
    3355          49 :     for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
    3356             :   }
    3357          14 :   return gerepilecopy(av,V);
    3358             : }

Generated by: LCOV version 1.13