Code coverage tests

This page documents the degree to which the PARI/GP source code is tested by our public test suite, distributed with the source distribution in directory src/test/. This is measured by the gcov utility; we then process gcov output using the lcov frond-end.

We test a few variants depending on Configure flags on the pari.math.u-bordeaux.fr machine (x86_64 architecture), and agregate them in the final report:

The target is 90% coverage for all mathematical modules (given that branches depending on DEBUGLEVEL or DEBUGMEM are not covered). This script is run to produce the results below.

LCOV - code coverage report
Current view: top level - basemath - base3.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.8.0 lcov report (development 19199-979caa2) Lines: 1369 1456 94.0 %
Date: 2016-07-24 07:10:30 Functions: 132 138 95.7 %
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     9927581 : get_tab(GEN nf, long *N)
      30             : {
      31     9927581 :   GEN tab = (typ(nf) == t_MAT)? nf: gel(nf,9);
      32     9927581 :   *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   429431303 : _mulii(GEN x, GEN y) {
      38  1084936103 :   return is_pm1(x)? (signe(x) < 0)? negi(y): y
      39   655504800 :                   : mulii(x, y);
      40             : }
      41             : 
      42             : GEN
      43        2310 : tablemul_ei_ej(GEN M, long i, long j)
      44             : {
      45             :   long N;
      46        2310 :   GEN tab = get_tab(M, &N);
      47        2310 :   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       16639 : tablemul_ei(GEN M, GEN x, long i)
      55             : {
      56             :   long j, k, N;
      57             :   GEN v, tab;
      58             : 
      59       16639 :   if (i==1) return gcopy(x);
      60       14063 :   tab = get_tab(M, &N);
      61       14063 :   if (typ(x) != t_COL) { v = zerocol(N); gel(v,i) = gcopy(x); return v; }
      62       14063 :   tab += (i-1)*N; v = cgetg(N+1,t_COL);
      63             :   /* wi . x = [ sum_j tab[k,j] x[j] ]_k */
      64      102431 :   for (k=1; k<=N; k++)
      65             :   {
      66       88368 :     pari_sp av = avma;
      67       88368 :     GEN s = gen_0;
      68      686490 :     for (j=1; j<=N; j++)
      69             :     {
      70      598122 :       GEN c = gcoeff(tab,k,j);
      71      598122 :       if (!gequal0(c)) s = gadd(s, gmul(c, gel(x,j)));
      72             :     }
      73       88368 :     gel(v,k) = gerepileupto(av,s);
      74             :   }
      75       14063 :   return v;
      76             : }
      77             : /* as tablemul_ei, assume x a ZV of correct length */
      78             : GEN
      79     8756865 : zk_ei_mul(GEN nf, GEN x, long i)
      80             : {
      81             :   long j, k, N;
      82             :   GEN v, tab;
      83             : 
      84     8756865 :   if (i==1) return ZC_copy(x);
      85     8756851 :   tab = get_tab(nf, &N); tab += (i-1)*N;
      86     8756851 :   v = cgetg(N+1,t_COL);
      87    65623121 :   for (k=1; k<=N; k++)
      88             :   {
      89    56866270 :     pari_sp av = avma;
      90    56866270 :     GEN s = gen_0;
      91   723852204 :     for (j=1; j<=N; j++)
      92             :     {
      93   666985934 :       GEN c = gcoeff(tab,k,j);
      94   666985934 :       if (signe(c)) s = addii(s, _mulii(c, gel(x,j)));
      95             :     }
      96    56866270 :     gel(v,k) = gerepileuptoint(av, s);
      97             :   }
      98     8756851 :   return v;
      99             : }
     100             : 
     101             : /* table of multiplication by wi in R[w1,..., wN] */
     102             : GEN
     103        4375 : ei_multable(GEN TAB, long i)
     104             : {
     105             :   long k,N;
     106        4375 :   GEN m, tab = get_tab(TAB, &N);
     107        4375 :   tab += (i-1)*N;
     108        4375 :   m = cgetg(N+1,t_MAT);
     109        4375 :   for (k=1; k<=N; k++) gel(m,k) = gel(tab,k);
     110        4375 :   return m;
     111             : }
     112             : 
     113             : GEN
     114     3562507 : zk_multable(GEN nf, GEN x)
     115             : {
     116     3562507 :   long i, l = lg(x);
     117     3562507 :   GEN mul = cgetg(l,t_MAT);
     118     3562507 :   gel(mul,1) = x; /* assume w_1 = 1 */
     119     3562507 :   for (i=2; i<l; i++) gel(mul,i) = zk_ei_mul(nf,x,i);
     120     3562507 :   return mul;
     121             : }
     122             : GEN
     123         707 : multable(GEN M, GEN x)
     124             : {
     125             :   long i, N;
     126             :   GEN mul;
     127         707 :   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     6860200 : zk_scalar_or_multable(GEN nf, GEN x)
     140             : {
     141     6860200 :   long tx = typ(x);
     142     6860200 :   if (tx == t_MAT || tx == t_INT) return x;
     143     3499535 :   x = nf_to_scalar_or_basis(nf, x);
     144     3499528 :   return (typ(x) == t_COL)? zk_multable(nf, x): x;
     145             : }
     146             : 
     147             : GEN
     148          42 : nftrace(GEN nf, GEN x)
     149             : {
     150          42 :   pari_sp av = avma;
     151          42 :   nf = checknf(nf);
     152          42 :   x = nf_to_scalar_or_basis(nf, x);
     153         105 :   x = (typ(x) == t_COL)? RgV_dotproduct(x, gel(nf_get_Tr(nf),1))
     154          63 :                        : gmulgs(x, nf_get_degree(nf));
     155          42 :   return gerepileupto(av, x);
     156             : }
     157             : GEN
     158         567 : rnfelttrace(GEN rnf, GEN x)
     159             : {
     160         567 :   pari_sp av = avma;
     161         567 :   checkrnf(rnf);
     162         567 :   x = rnfeltabstorel(rnf, x);
     163        1344 :   x = (typ(x) == t_POLMOD)? rnfeltdown(rnf, gtrace(x))
     164         868 :                           : gmulgs(x, rnf_get_degree(rnf));
     165         476 :   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       18851 : nfnorm(GEN nf, GEN x)
     181             : {
     182       18851 :   pari_sp av = avma;
     183       18851 :   nf = checknf(nf);
     184       18851 :   if (typ(x) == t_MAT) return famat_norm(nf, x);
     185       18844 :   x = nf_to_scalar_or_alg(nf, x);
     186       54187 :   x = (typ(x) == t_POL)? RgXQ_norm(x, nf_get_pol(nf))
     187       35343 :                        : gpowgs(x, nf_get_degree(nf));
     188       18844 :   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       24983 : nfadd(GEN nf, GEN x, GEN y)
     205             : {
     206       24983 :   pari_sp av = avma;
     207             :   GEN z;
     208             : 
     209       24983 :   nf = checknf(nf);
     210       24983 :   x = nf_to_scalar_or_basis(nf, x);
     211       24983 :   y = nf_to_scalar_or_basis(nf, y);
     212       24983 :   if (typ(x) != t_COL) {
     213        7070 :     if (typ(y) == t_COL) z = RgC_Rg_add(y, x);
     214             :     else {
     215        3143 :       long N = nf_get_degree(nf);
     216        3143 :       z = zerocol(N); gel(z,1) = gadd(x,y);
     217             :     }
     218             :   } else {
     219       17913 :     if (typ(y) != t_COL) z = RgC_Rg_add(x, y);
     220        7091 :     else z = RgC_add(x, y);
     221             :   }
     222       24983 :   return gerepileupto(av, z);
     223             : }
     224             : /* x - y in nf */
     225             : GEN
     226           7 : nfsub(GEN nf, GEN x, GEN y)
     227             : {
     228           7 :   pari_sp av = avma;
     229             :   GEN z;
     230             : 
     231           7 :   nf = checknf(nf);
     232           7 :   x = nf_to_scalar_or_basis(nf, x);
     233           7 :   y = nf_to_scalar_or_basis(nf, y);
     234           7 :   if (typ(x) != t_COL) {
     235           0 :     if (typ(y) == t_COL) z = RgC_Rg_sub(y, x);
     236             :     else {
     237           0 :       long N = nf_get_degree(nf);
     238           0 :       z = zerocol(N); gel(z,1) = gsub(x,y);
     239             :     }
     240             :   } else {
     241           7 :     if (typ(y) != t_COL) z = RgC_Rg_sub(x, y);
     242           7 :     else z = RgC_sub(x, y);
     243             :   }
     244           7 :   return gerepileupto(av, z);
     245             : }
     246             : 
     247             : /* product of x and y in nf */
     248             : GEN
     249      577560 : nfmul(GEN nf, GEN x, GEN y)
     250             : {
     251             :   GEN z;
     252      577560 :   pari_sp av = avma;
     253             : 
     254      577560 :   if (x == y) return nfsqr(nf,x);
     255             : 
     256      577490 :   nf = checknf(nf);
     257      577490 :   x = nf_to_scalar_or_basis(nf, x);
     258      577490 :   y = nf_to_scalar_or_basis(nf, y);
     259      577490 :   if (typ(x) != t_COL)
     260             :   {
     261      445522 :     if (typ(y) == t_COL) z = RgC_Rg_mul(y, x);
     262             :     else {
     263      394520 :       long N = nf_get_degree(nf);
     264      394520 :       z = zerocol(N); gel(z,1) = gmul(x,y);
     265             :     }
     266             :   }
     267             :   else
     268             :   {
     269      131968 :     if (typ(y) != t_COL) z = RgC_Rg_mul(x, y);
     270             :     else {
     271             :       GEN dx, dy;
     272       66784 :       x = Q_remove_denom(x, &dx);
     273       66784 :       y = Q_remove_denom(y, &dy);
     274       66784 :       z = nfmuli(nf,x,y);
     275       66784 :       dx = mul_denom(dx,dy);
     276       66784 :       if (dx) z = RgC_Rg_div(z, dx);
     277             :     }
     278             :   }
     279      577490 :   return gerepileupto(av, z);
     280             : }
     281             : /* square of x in nf */
     282             : GEN
     283       11641 : nfsqr(GEN nf, GEN x)
     284             : {
     285       11641 :   pari_sp av = avma;
     286             :   GEN z;
     287             : 
     288       11641 :   nf = checknf(nf);
     289       11641 :   x = nf_to_scalar_or_basis(nf, x);
     290       11641 :   if (typ(x) != t_COL)
     291             :   {
     292       10073 :     long N = nf_get_degree(nf);
     293       10073 :     z = zerocol(N); gel(z,1) = gsqr(x);
     294             :   }
     295             :   else
     296             :   {
     297             :     GEN dx;
     298        1568 :     x = Q_remove_denom(x, &dx);
     299        1568 :     z = nfsqri(nf,x);
     300        1568 :     if (dx) z = RgC_Rg_div(z, sqri(dx));
     301             :   }
     302       11641 :   return gerepileupto(av, z);
     303             : }
     304             : 
     305             : GEN
     306      189317 : nfC_nf_mul(GEN nf, GEN v, GEN x)
     307             : {
     308             :   long tx, l, i;
     309      189317 :   GEN y, dx = NULL;
     310             : 
     311      189317 :   x = nf_to_scalar_or_basis(nf, x);
     312      189317 :   tx = typ(x);
     313      189317 :   if (tx != t_COL) {
     314       53230 :     if (tx == t_INT)
     315             :     {
     316       51928 :       long s = signe(x);
     317       51928 :       if (!s) return zerocol(lg(v)-1);
     318       49128 :       if (is_pm1(x)) return s > 0? leafcopy(v): RgC_neg(v);
     319             :     }
     320       14667 :     l = lg(v); y = cgetg(l, t_COL);
     321      112817 :     for (i=1; i < l; i++)
     322             :     {
     323       98150 :       GEN c = gel(v,i);
     324       98150 :       if (typ(c) != t_COL) c = gmul(c, x); else c = RgC_Rg_mul(c, x);
     325       98150 :       gel(y,i) = c;
     326             :     }
     327             :   }
     328             :   else
     329             :   {
     330      136087 :     x = Q_remove_denom(x, &dx);
     331      136087 :     x = zk_multable(nf, x);
     332      136087 :     l = lg(v); y = cgetg(l, t_COL);
     333      539032 :     for (i=1; i < l; i++)
     334             :     {
     335      402945 :       GEN c = gel(v,i);
     336      402945 :       if (typ(c)!=t_COL) {
     337       85341 :         if (!isintzero(c)) c = RgC_Rg_mul(gel(x,1), c);
     338             :       } else {
     339      317604 :         c = RgM_RgC_mul(x,c);
     340      317604 :         if (QV_isscalar(c)) c = gel(c,1);
     341             :       }
     342      402945 :       gel(y,i) = c;
     343             :     }
     344             :   }
     345      150754 :   return dx? RgC_Rg_div(y, dx): y;
     346             : }
     347             : static GEN
     348        3318 : mulbytab(GEN M, GEN c)
     349        3318 : { return typ(c) == t_COL? RgM_RgC_mul(M,c): RgC_Rg_mul(gel(M,1), c); }
     350             : GEN
     351         707 : tablemulvec(GEN M, GEN x, GEN v)
     352             : {
     353             :   long l, i;
     354             :   GEN y;
     355             : 
     356         707 :   if (typ(x) == t_COL && RgV_isscalar(x))
     357             :   {
     358           0 :     x = gel(x,1);
     359           0 :     return typ(v) == t_POL? RgX_Rg_mul(v,x): RgV_Rg_mul(v,x);
     360             :   }
     361         707 :   x = multable(M, x); /* multiplication table by x */
     362         707 :   y = cgetg_copy(v, &l);
     363         707 :   if (typ(v) == t_POL)
     364             :   {
     365         707 :     y[1] = v[1];
     366         707 :     for (i=2; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     367         707 :     y = normalizepol(y);
     368             :   }
     369             :   else
     370             :   {
     371           0 :     for (i=1; i < l; i++) gel(y,i) = mulbytab(x, gel(v,i));
     372             :   }
     373         707 :   return y;
     374             : }
     375             : 
     376             : /* inverse of x in nf */
     377             : GEN
     378      152712 : nfinv(GEN nf, GEN x)
     379             : {
     380      152712 :   pari_sp av = avma;
     381             :   GEN T, z;
     382             : 
     383      152712 :   nf = checknf(nf); T = nf_get_pol(nf);
     384      152712 :   x = nf_to_scalar_or_alg(nf, x);
     385      152712 :   if (typ(x) == t_POL)
     386        4613 :     z = poltobasis(nf, QXQ_inv(x, T));
     387             :   else {
     388      148099 :     z = zerocol(degpol(T)); gel(z,1) = ginv(x);
     389             :   }
     390      152712 :   return gerepileupto(av, z);
     391             : }
     392             : 
     393             : /* quotient of x and y in nf */
     394             : GEN
     395       12046 : nfdiv(GEN nf, GEN x, GEN y)
     396             : {
     397       12046 :   pari_sp av = avma;
     398             :   GEN T, z;
     399             : 
     400       12046 :   nf = checknf(nf); T = nf_get_pol(nf);
     401       12046 :   y = nf_to_scalar_or_alg(nf, y);
     402       11780 :   if (typ(y) != t_POL) {
     403        2415 :     x = nf_to_scalar_or_basis(nf, x);
     404        2415 :     if (typ(x) == t_COL) z = RgC_Rg_div(x, y);
     405             :     else {
     406        1036 :       z = zerocol(degpol(T)); gel(z,1) = gdiv(x,y);
     407             :     }
     408             :   }
     409             :   else
     410             :   {
     411        9365 :     x = nf_to_scalar_or_alg(nf, x);
     412        9323 :     z = QXQ_inv(y, T);
     413        9323 :     z = (typ(x) == t_POL)? RgXQ_mul(z, x, T): RgX_Rg_mul(z, x);
     414        9323 :     z = poltobasis(nf, z);
     415             :   }
     416       11710 :   return gerepileupto(av, z);
     417             : }
     418             : 
     419             : /* product of INTEGERS (t_INT or ZC) x and y in nf
     420             :  * compute xy as ( sum_i x_i sum_j y_j m^{i,j}_k )_k */
     421             : GEN
     422      464153 : nfmuli(GEN nf, GEN x, GEN y)
     423             : {
     424             :   long i, j, k, N;
     425      464153 :   GEN s, v, TAB = get_tab(nf, &N);
     426             : 
     427      464153 :   if (typ(x) == t_INT)
     428             :   {
     429       67651 :     if (typ(y) == t_INT) return scalarcol(mulii(x,y), N);
     430        5635 :     return ZC_Z_mul(y, x);
     431             :   }
     432      396502 :   if (typ(y) == t_INT) return ZC_Z_mul(x, y);
     433             :   /* both x and y are ZV */
     434      372349 :   v = cgetg(N+1,t_COL);
     435     2329853 :   for (k=1; k<=N; k++)
     436             :   {
     437     1957504 :     pari_sp av = avma;
     438     1957504 :     GEN TABi = TAB;
     439     1957504 :     if (k == 1)
     440      372349 :       s = mulii(gel(x,1),gel(y,1));
     441             :     else
     442     3170310 :       s = addii(mulii(gel(x,1),gel(y,k)),
     443     3170310 :                 mulii(gel(x,k),gel(y,1)));
     444    15285654 :     for (i=2; i<=N; i++)
     445             :     {
     446    13328150 :       GEN t, xi = gel(x,i);
     447    13328150 :       TABi += N;
     448    13328150 :       if (!signe(xi)) continue;
     449             : 
     450     8562349 :       t = NULL;
     451   115054519 :       for (j=2; j<=N; j++)
     452             :       {
     453   106492170 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     454   106492170 :         if (!signe(c)) continue;
     455    44748339 :         p1 = _mulii(c, gel(y,j));
     456    44748339 :         t = t? addii(t, p1): p1;
     457             :       }
     458     8562349 :       if (t) s = addii(s, mulii(xi, t));
     459             :     }
     460     1957504 :     gel(v,k) = gerepileuptoint(av,s);
     461             :   }
     462      372349 :   return v;
     463             : }
     464             : /* square of INTEGER (t_INT or ZC) x in nf */
     465             : GEN
     466      685829 : nfsqri(GEN nf, GEN x)
     467             : {
     468             :   long i, j, k, N;
     469      685829 :   GEN s, v, TAB = get_tab(nf, &N);
     470             : 
     471      685829 :   if (typ(x) == t_INT) return scalarcol(sqri(x), N);
     472      685829 :   v = cgetg(N+1,t_COL);
     473     5622900 :   for (k=1; k<=N; k++)
     474             :   {
     475     4937071 :     pari_sp av = avma;
     476     4937071 :     GEN TABi = TAB;
     477     4937071 :     if (k == 1)
     478      685829 :       s = sqri(gel(x,1));
     479             :     else
     480     4251242 :       s = shifti(mulii(gel(x,1),gel(x,k)), 1);
     481    58091791 :     for (i=2; i<=N; i++)
     482             :     {
     483    53154720 :       GEN p1, c, t, xi = gel(x,i);
     484    53154720 :       TABi += N;
     485    53154720 :       if (!signe(xi)) continue;
     486             : 
     487    18369618 :       c = gcoeff(TABi, k, i);
     488    18369618 :       t = signe(c)? _mulii(c,xi): NULL;
     489   261100649 :       for (j=i+1; j<=N; j++)
     490             :       {
     491   242731031 :         c = gcoeff(TABi, k, j);
     492   242731031 :         if (!signe(c)) continue;
     493   127981637 :         p1 = _mulii(c, shifti(gel(x,j),1));
     494   127981637 :         t = t? addii(t, p1): p1;
     495             :       }
     496    18369618 :       if (t) s = addii(s, mulii(xi, t));
     497             :     }
     498     4937071 :     gel(v,k) = gerepileuptoint(av,s);
     499             :   }
     500      685829 :   return v;
     501             : }
     502             : 
     503             : /* both x and y are RgV */
     504             : GEN
     505           0 : tablemul(GEN TAB, GEN x, GEN y)
     506             : {
     507             :   long i, j, k, N;
     508             :   GEN s, v;
     509           0 :   if (typ(x) != t_COL) return gmul(x, y);
     510           0 :   if (typ(y) != t_COL) return gmul(y, x);
     511           0 :   N = lg(x)-1;
     512           0 :   v = cgetg(N+1,t_COL);
     513           0 :   for (k=1; k<=N; k++)
     514             :   {
     515           0 :     pari_sp av = avma;
     516           0 :     GEN TABi = TAB;
     517           0 :     if (k == 1)
     518           0 :       s = gmul(gel(x,1),gel(y,1));
     519             :     else
     520           0 :       s = gadd(gmul(gel(x,1),gel(y,k)),
     521           0 :                gmul(gel(x,k),gel(y,1)));
     522           0 :     for (i=2; i<=N; i++)
     523             :     {
     524           0 :       GEN t, xi = gel(x,i);
     525           0 :       TABi += N;
     526           0 :       if (gequal0(xi)) continue;
     527             : 
     528           0 :       t = NULL;
     529           0 :       for (j=2; j<=N; j++)
     530             :       {
     531           0 :         GEN p1, c = gcoeff(TABi, k, j); /* m^{i,j}_k */
     532           0 :         if (gequal0(c)) continue;
     533           0 :         p1 = gmul(c, gel(y,j));
     534           0 :         t = t? gadd(t, p1): p1;
     535             :       }
     536           0 :       if (t) s = gadd(s, gmul(xi, t));
     537             :     }
     538           0 :     gel(v,k) = gerepileupto(av,s);
     539             :   }
     540           0 :   return v;
     541             : }
     542             : GEN
     543        5726 : tablesqr(GEN TAB, GEN x)
     544             : {
     545             :   long i, j, k, N;
     546             :   GEN s, v;
     547             : 
     548        5726 :   if (typ(x) != t_COL) return gsqr(x);
     549        5726 :   N = lg(x)-1;
     550        5726 :   v = cgetg(N+1,t_COL);
     551             : 
     552       42784 :   for (k=1; k<=N; k++)
     553             :   {
     554       37058 :     pari_sp av = avma;
     555       37058 :     GEN TABi = TAB;
     556       37058 :     if (k == 1)
     557        5726 :       s = gsqr(gel(x,1));
     558             :     else
     559       31332 :       s = gmul2n(gmul(gel(x,1),gel(x,k)), 1);
     560      253372 :     for (i=2; i<=N; i++)
     561             :     {
     562      216314 :       GEN p1, c, t, xi = gel(x,i);
     563      216314 :       TABi += N;
     564      216314 :       if (gequal0(xi)) continue;
     565             : 
     566       71015 :       c = gcoeff(TABi, k, i);
     567       71015 :       t = !gequal0(c)? gmul(c,xi): NULL;
     568      318045 :       for (j=i+1; j<=N; j++)
     569             :       {
     570      247030 :         c = gcoeff(TABi, k, j);
     571      247030 :         if (gequal0(c)) continue;
     572      131278 :         p1 = gmul(gmul2n(c,1), gel(x,j));
     573      131278 :         t = t? gadd(t, p1): p1;
     574             :       }
     575       71015 :       if (t) s = gadd(s, gmul(xi, t));
     576             :     }
     577       37058 :     gel(v,k) = gerepileupto(av,s);
     578             :   }
     579        5726 :   return v;
     580             : }
     581             : 
     582             : static GEN
     583       45331 : _mul(void *data, GEN x, GEN y) { return nfmuli((GEN)data,x,y); }
     584             : static GEN
     585      264778 : _sqr(void *data, GEN x) { return nfsqri((GEN)data,x); }
     586             : 
     587             : /* Compute z^n in nf, left-shift binary powering */
     588             : GEN
     589      205435 : nfpow(GEN nf, GEN z, GEN n)
     590             : {
     591      205435 :   pari_sp av = avma;
     592             :   long s, N;
     593             :   GEN x, cx, T;
     594             : 
     595      205435 :   if (typ(n)!=t_INT) pari_err_TYPE("nfpow",n);
     596      205435 :   nf = checknf(nf); T = nf_get_pol(nf); N = degpol(T);
     597      205435 :   s = signe(n); if (!s) return scalarcol_shallow(gen_1,N);
     598      205435 :   x = nf_to_scalar_or_basis(nf, z);
     599      205435 :   if (typ(x) != t_COL) { GEN y = zerocol(N); gel(y,1) = powgi(x,n); return y; }
     600      204245 :   if (s < 0) { /* simplified nfinv */
     601        1876 :     x = nf_to_scalar_or_alg(nf, z);
     602        1876 :     x = poltobasis(nf, QXQ_inv(x, T));
     603        1876 :     n = absi(n);
     604             :   }
     605      204245 :   x = primitive_part(x, &cx);
     606      204245 :   x = gen_pow(x, n, (void*)nf, _sqr, _mul);
     607      204245 :   if (cx) x = RgC_Rg_mul(x, powgi(cx, n));
     608      204245 :   return av==avma? gcopy(x): gerepileupto(av,x);
     609             : }
     610             : /* Compute z^n in nf, left-shift binary powering */
     611             : GEN
     612        4767 : nfpow_u(GEN nf, GEN z, ulong n)
     613             : {
     614        4767 :   pari_sp av = avma;
     615             :   long N;
     616             :   GEN x, cx, T;
     617             : 
     618        4767 :   nf = checknf(nf); T = nf_get_pol(nf); N = degpol(T);
     619        4767 :   if (!n) return scalarcol_shallow(gen_1,N);
     620        4767 :   x = nf_to_scalar_or_basis(nf, z);
     621        4767 :   if (typ(x) != t_COL) { GEN y = zerocol(N); gel(y,1) = gpowgs(x,n); return y; }
     622        4025 :   x = primitive_part(x, &cx);
     623        4025 :   x = gen_powu(x, n, (void*)nf, _sqr, _mul);
     624        4025 :   if (cx) x = RgC_Rg_mul(x, powgi(cx, utoipos(n)));
     625        4025 :   return av==avma? gcopy(x): gerepileupto(av,x);
     626             : }
     627             : 
     628             : typedef struct {
     629             :   GEN nf, p;
     630             :   long I;
     631             : } eltmod_muldata;
     632             : 
     633             : static GEN
     634      136746 : sqr_mod(void *data, GEN x)
     635             : {
     636      136746 :   eltmod_muldata *D = (eltmod_muldata*)data;
     637      136746 :   return FpC_red(nfsqri(D->nf, x), D->p);
     638             : }
     639             : static GEN
     640       64044 : ei_msqr_mod(void *data, GEN x)
     641             : {
     642       64044 :   GEN x2 = sqr_mod(data, x);
     643       64044 :   eltmod_muldata *D = (eltmod_muldata*)data;
     644       64044 :   return FpC_red(zk_ei_mul(D->nf, x2, D->I), D->p);
     645             : }
     646             : 
     647             : /* x = I-th vector of the Z-basis of Z_K, in Z^n, compute lift(x^n mod p) */
     648             : GEN
     649       91531 : pow_ei_mod_p(GEN nf, long I, GEN n, GEN p)
     650             : {
     651       91531 :   pari_sp av = avma;
     652             :   eltmod_muldata D;
     653             :   long s,N;
     654             :   GEN y;
     655             : 
     656       91531 :   if (typ(n) != t_INT) pari_err_TYPE("nfpow",n);
     657       91531 :   nf = checknf(nf); N = nf_get_degree(nf);
     658       91531 :   s = signe(n);
     659       91531 :   if (s < 0) pari_err_IMPL("negative power in pow_ei_mod_p");
     660       91531 :   if (!s || I == 1) return scalarcol_shallow(gen_1,N);
     661       76787 :   D.nf = nf;
     662       76787 :   D.p = p;
     663       76787 :   D.I = I;
     664       76787 :   y = gen_pow_fold(col_ei(N, I), n, (void*)&D, &sqr_mod, &ei_msqr_mod);
     665       76787 :   return gerepileupto(av,y);
     666             : }
     667             : 
     668             : /* valuation of integral x (ZV), with resp. to prime ideal pr */
     669             : long
     670     3355381 : ZC_nfvalrem(GEN nf, GEN x, GEN pr, GEN *newx)
     671             : {
     672             :   long i, v, l;
     673     3355381 :   GEN r, y, p = pr_get_p(pr), mul = zk_scalar_or_multable(nf, pr_get_tau(pr));
     674             : 
     675             :   /* p inert */
     676     3355381 :   if (typ(mul) == t_INT) return newx? ZV_pvalrem(x, p, newx):ZV_pval(x, p);
     677     3272844 :   y = cgetg_copy(x, &l); /* will hold the new x */
     678     3272844 :   x = leafcopy(x);
     679     4850215 :   for(v=0;; v++)
     680             :   {
     681    16559631 :     for (i=1; i<l; i++)
     682             :     { /* is (x.b)[i] divisible by p ? */
     683    14982260 :       gel(y,i) = dvmdii(ZMrow_ZC_mul(mul,x,i),p,&r);
     684    14982260 :       if (r != gen_0) { if (newx) *newx = x; return v; }
     685             :     }
     686     1577371 :     swap(x, y);
     687     1577371 :   }
     688             : }
     689             : long
     690     3150070 : ZC_nfval(GEN nf, GEN x, GEN P)
     691     3150070 : { return ZC_nfvalrem(nf, x, P, NULL); }
     692             : 
     693             : /* v_P(x) != 0, x a ZV. Simpler version of ZC_nfvalrem */
     694             : int
     695        8211 : ZC_prdvd(GEN nf, GEN x, GEN P)
     696             : {
     697        8211 :   pari_sp av = avma;
     698             :   long i, l;
     699        8211 :   GEN p = pr_get_p(P), mul = zk_scalar_or_multable(nf, pr_get_tau(P));
     700        8211 :   if (typ(mul) == t_INT) return ZV_Z_dvd(x, p);
     701        8162 :   l = lg(x);
     702       44205 :   for (i=1; i<l; i++)
     703       40698 :     if (remii(ZMrow_ZC_mul(mul,x,i), p) != gen_0) { avma = av; return 0; }
     704        3507 :   avma = av; return 1;
     705             : }
     706             : 
     707             : int
     708          28 : pr_equal(GEN nf, GEN P, GEN Q)
     709             : {
     710          28 :   GEN gQ, p = pr_get_p(P);
     711          28 :   long e = pr_get_e(P), f = pr_get_f(P), n;
     712          28 :   if (!equalii(p, pr_get_p(Q)) || e != pr_get_e(Q) || f != pr_get_f(Q))
     713          14 :     return 0;
     714          14 :   gQ = pr_get_gen(Q); n = lg(gQ)-1;
     715          14 :   if (2*e*f > n) return 1; /* room for only one such pr */
     716           7 :   return ZV_equal(pr_get_gen(P), gQ) || ZC_prdvd(nf, gQ, P);
     717             : }
     718             : 
     719             : long
     720      919590 : nfval(GEN nf, GEN x, GEN pr)
     721             : {
     722      919590 :   pari_sp av = avma;
     723             :   long w, e;
     724             :   GEN cx, p;
     725             : 
     726      919590 :   if (gequal0(x)) return LONG_MAX;
     727      913514 :   nf = checknf(nf);
     728      913514 :   checkprid(pr);
     729      913514 :   p = pr_get_p(pr);
     730      913514 :   e = pr_get_e(pr);
     731      913514 :   x = nf_to_scalar_or_basis(nf, x);
     732      913514 :   if (typ(x) != t_COL) return e*Q_pval(x,p);
     733      398517 :   x = Q_primitive_part(x, &cx);
     734      398517 :   w = ZC_nfval(nf,x,pr);
     735      398517 :   if (cx) w += e*Q_pval(cx,p);
     736      398517 :   avma = av; return w;
     737             : }
     738             : 
     739             : /* want to write p^v = uniformizer^(e*v) * z^v, z coprime to pr */
     740             : /* z := tau^e / p^(e-1), algebraic integer coprime to pr; return z^v */
     741             : static GEN
     742        3808 : powp(GEN nf, GEN pr, long v)
     743             : {
     744             :   GEN b, z;
     745             :   long e;
     746        3808 :   if (!v) return gen_1;
     747        3794 :   b = pr_get_tau(pr);
     748        3794 :   if (typ(b) == t_INT) return gen_1;
     749         574 :   e = pr_get_e(pr);
     750         574 :   z = gel(b,1);
     751         574 :   if (e != 1) z = gdiv(nfpow_u(nf, z, e), powiu(pr_get_p(pr),e-1));
     752         574 :   return nfpow_u(nf, z, v);
     753             : }
     754             : long
     755        9639 : nfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     756             : {
     757        9639 :   pari_sp av = avma;
     758             :   long w, e;
     759             :   GEN cx, p, t;
     760             : 
     761        9639 :   if (!py) return nfval(nf,x,pr);
     762        9520 :   if (gequal0(x)) { *py = gcopy(x); return LONG_MAX; }
     763        9506 :   nf = checknf(nf);
     764        9506 :   checkprid(pr);
     765        9506 :   p = pr_get_p(pr);
     766        9506 :   e = pr_get_e(pr);
     767        9506 :   x = nf_to_scalar_or_basis(nf, x);
     768        9506 :   if (typ(x) != t_COL) {
     769        3367 :     w = Q_pvalrem(x,p, py);
     770        3367 :     if (!w) { *py = gerepilecopy(av, x); return 0; }
     771        3290 :     *py = gerepileupto(av, gmul(powp(nf, pr, w), *py));
     772        3290 :     return e*w;
     773             :   }
     774        6139 :   x = Q_primitive_part(x, &cx);
     775        6139 :   w = ZC_nfvalrem(nf,x,pr, py);
     776        6139 :   if (cx)
     777             :   {
     778         518 :     long v = Q_pvalrem(cx,p, &t);
     779         518 :     *py = nfmul(nf, *py, gmul(powp(nf,pr,v), t));
     780         518 :     *py = gerepileupto(av, *py);
     781         518 :     w += e*v;
     782             :   }
     783             :   else
     784        5621 :     *py = gerepilecopy(av, *py);
     785        6139 :   return w;
     786             : }
     787             : GEN
     788         147 : gpnfvalrem(GEN nf, GEN x, GEN pr, GEN *py)
     789             : {
     790         147 :   long v = nfvalrem(nf,x,pr,py);
     791         147 :   return v == LONG_MAX? mkoo(): stoi(v);
     792             : }
     793             : 
     794             : GEN
     795      150332 : coltoalg(GEN nf, GEN x)
     796             : {
     797      150332 :   return mkpolmod( coltoliftalg(nf, x), nf_get_pol(nf) );
     798             : }
     799             : 
     800             : GEN
     801      609637 : basistoalg(GEN nf, GEN x)
     802             : {
     803             :   GEN z, T;
     804             : 
     805      609637 :   nf = checknf(nf);
     806      609637 :   switch(typ(x))
     807             :   {
     808             :     case t_COL: {
     809      143864 :       pari_sp av = avma;
     810      143864 :       return gerepilecopy(av, coltoalg(nf, x));
     811             :     }
     812             :     case t_POLMOD:
     813         119 :       T = nf_get_pol(nf);
     814         119 :       if (!RgX_equal_var(T,gel(x,1)))
     815           0 :         pari_err_MODULUS("basistoalg", T,gel(x,1));
     816         119 :       return gcopy(x);
     817             :     case t_POL:
     818         511 :       T = nf_get_pol(nf);
     819         511 :       if (varn(T) != varn(x)) pari_err_VAR("basistoalg",x,T);
     820         511 :       z = cgetg(3,t_POLMOD);
     821         511 :       gel(z,1) = ZX_copy(T);
     822         511 :       gel(z,2) = RgX_rem(x, T); return z;
     823             :     case t_INT:
     824             :     case t_FRAC:
     825      465143 :       T = nf_get_pol(nf);
     826      465143 :       z = cgetg(3,t_POLMOD);
     827      465143 :       gel(z,1) = ZX_copy(T);
     828      465143 :       gel(z,2) = gcopy(x); return z;
     829             :     default:
     830           0 :       pari_err_TYPE("basistoalg",x);
     831           0 :       return NULL; /* not reached */
     832             :   }
     833             : }
     834             : 
     835             : /* Assume nf is a genuine nf. */
     836             : GEN
     837     8622291 : nf_to_scalar_or_basis(GEN nf, GEN x)
     838             : {
     839     8622291 :   switch(typ(x))
     840             :   {
     841             :     case t_INT: case t_FRAC:
     842     1404381 :       return x;
     843             :     case t_POLMOD:
     844     1106812 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_basis");
     845     1106742 :       if (typ(x) != t_POL) return x;
     846             :       /* fall through */
     847             :     case t_POL:
     848             :     {
     849      553272 :       GEN T = nf_get_pol(nf);
     850      553272 :       long l = lg(x);
     851      553272 :       if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_basis", x,T);
     852      553216 :       if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
     853      553216 :       if (l == 2) return gen_0;
     854      528261 :       if (l == 3) return gel(x,2);
     855      282911 :       return poltobasis(nf,x);
     856             :     }
     857             :     case t_COL:
     858     6071885 :       if (lg(x) != lg(nf_get_zk(nf))) break;
     859     6071815 :       return QV_isscalar(x)? gel(x,1): x;
     860             :   }
     861          70 :   pari_err_TYPE("nf_to_scalar_or_basis",x);
     862           0 :   return NULL; /* not reached */
     863             : }
     864             : /* Let x be a polynomial with coefficients in Q or nf. Return the same
     865             :  * polynomial with coefficients expressed as vectors (on the integral basis).
     866             :  * No consistency checks, not memory-clean. */
     867             : GEN
     868        2555 : RgX_to_nfX(GEN nf, GEN x)
     869             : {
     870             :   long i, l;
     871        2555 :   GEN y = cgetg_copy(x, &l); y[1] = x[1];
     872        2555 :   for (i=2; i<l; i++) gel(y,i) = nf_to_scalar_or_basis(nf, gel(x,i));
     873        2555 :   return y;
     874             : }
     875             : 
     876             : /* Assume nf is a genuine nf. */
     877             : GEN
     878      218410 : nf_to_scalar_or_alg(GEN nf, GEN x)
     879             : {
     880      218410 :   switch(typ(x))
     881             :   {
     882             :     case t_INT: case t_FRAC:
     883       56531 :       return x;
     884             :     case t_POLMOD:
     885       99582 :       x = checknfelt_mod(nf,x,"nf_to_scalar_or_alg");
     886       99582 :       if (typ(x) != t_POL) return x;
     887             :       /* fall through */
     888             :     case t_POL:
     889             :     {
     890       14168 :       GEN T = nf_get_pol(nf);
     891       14168 :       long l = lg(x);
     892       14168 :       if (varn(x) != varn(T)) pari_err_VAR("nf_to_scalar_or_alg", x,T);
     893       14028 :       if (l >= lg(T)) { x = RgX_rem(x, T); l = lg(x); }
     894       14028 :       if (l == 2) return gen_0;
     895       14028 :       if (l == 3) return gel(x,2);
     896       13426 :       return x;
     897             :     }
     898             :     case t_COL:
     899       49977 :       if (lg(x) != lg(nf_get_zk(nf))) break;
     900       49809 :       return QV_isscalar(x)? gel(x,1): coltoliftalg(nf, x);
     901             :   }
     902         217 :   pari_err_TYPE("nf_to_scalar_or_alg",x);
     903           0 :   return NULL; /* not reached */
     904             : }
     905             : 
     906             : /* gmul(A, RgX_to_RgC(x)), A t_MAT (or t_VEC) of compatible dimensions */
     907             : GEN
     908     1986241 : mulmat_pol(GEN A, GEN x)
     909             : {
     910             :   long i,l;
     911             :   GEN z;
     912     1986241 :   if (typ(x) != t_POL) return gmul(x,gel(A,1)); /* scalar */
     913     1986122 :   l=lg(x)-1; if (l == 1) return typ(A)==t_VEC? gen_0: zerocol(nbrows(A));
     914     1983938 :   x++; z = gmul(gel(x,1), gel(A,1));
     915     8129774 :   for (i=2; i<l ; i++)
     916     6145836 :     if (!gequal0(gel(x,i))) z = gadd(z, gmul(gel(x,i), gel(A,i)));
     917     1983938 :   return z;
     918             : }
     919             : 
     920             : /* x a t_POL, nf a genuine nf. No garbage collecting. No check.  */
     921             : GEN
     922     1829219 : poltobasis(GEN nf, GEN x)
     923             : {
     924     1829219 :   GEN P = nf_get_pol(nf);
     925     1829219 :   if (varn(x) != varn(P)) pari_err_VAR( "poltobasis", x,P);
     926     1829009 :   if (degpol(x) >= degpol(P)) x = RgX_rem(x,P);
     927     1829009 :   return mulmat_pol(nf_get_invzk(nf), x);
     928             : }
     929             : 
     930             : GEN
     931       64300 : algtobasis(GEN nf, GEN x)
     932             : {
     933             :   pari_sp av;
     934             : 
     935       64300 :   nf = checknf(nf);
     936       64300 :   switch(typ(x))
     937             :   {
     938             :     case t_POLMOD:
     939       39333 :       if (!RgX_equal_var(nf_get_pol(nf),gel(x,1)))
     940           7 :         pari_err_MODULUS("algtobasis", nf_get_pol(nf),gel(x,1));
     941       39326 :       x = gel(x,2);
     942       39326 :       switch(typ(x))
     943             :       {
     944             :         case t_INT:
     945        2877 :         case t_FRAC: return scalarcol(x, nf_get_degree(nf));
     946             :         case t_POL:
     947       36449 :           av = avma;
     948       36449 :           return gerepileupto(av,poltobasis(nf,x));
     949             :       }
     950           0 :       break;
     951             : 
     952             :     case t_POL:
     953        8328 :       av = avma;
     954        8328 :       return gerepileupto(av,poltobasis(nf,x));
     955             : 
     956             :     case t_COL:
     957        8225 :       if (lg(x)-1 != nf_get_degree(nf)) pari_err_DIM("nfalgtobasis");
     958        8225 :       return gcopy(x);
     959             : 
     960             :     case t_INT:
     961        8414 :     case t_FRAC: return scalarcol(x, nf_get_degree(nf));
     962             :   }
     963           0 :   pari_err_TYPE("algtobasis",x);
     964           0 :   return NULL; /* not reached */
     965             : }
     966             : 
     967             : GEN
     968       35637 : rnfbasistoalg(GEN rnf,GEN x)
     969             : {
     970       35637 :   const char *f = "rnfbasistoalg";
     971             :   long lx, i;
     972       35637 :   pari_sp av = avma;
     973             :   GEN z, nf, relpol, T;
     974             : 
     975       35637 :   checkrnf(rnf);
     976       35637 :   nf = rnf_get_nf(rnf);
     977       35637 :   T = nf_get_pol(nf);
     978       35637 :   relpol = QXQX_to_mod_shallow(rnf_get_pol(rnf), T);
     979       35637 :   switch(typ(x))
     980             :   {
     981             :     case t_COL:
     982         798 :       z = cgetg_copy(x, &lx);
     983        2338 :       for (i=1; i<lx; i++)
     984             :       {
     985        1589 :         GEN c = nf_to_scalar_or_alg(nf, gel(x,i));
     986        1540 :         if (typ(c) == t_POL) c = mkpolmod(c,T);
     987        1540 :         gel(z,i) = c;
     988             :       }
     989         749 :       z = RgV_RgC_mul(gel(rnf_get_zk(rnf),1), z);
     990         686 :       return gerepileupto(av, gmodulo(z,relpol));
     991             : 
     992             :     case t_POLMOD:
     993       23688 :       x = polmod_nffix(f, rnf, x, 0);
     994       23478 :       if (typ(x) != t_POL) break;
     995        9415 :       retmkpolmod(RgX_copy(x), RgX_copy(relpol));
     996             :     case t_POL:
     997         819 :       if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = gmodulo(x,T); break; }
     998         595 :       if (varn(x) == varn(relpol))
     999             :       {
    1000         546 :         x = RgX_nffix(f,nf_get_pol(nf),x,0);
    1001         546 :         return gmodulo(x, relpol);
    1002             :       }
    1003          49 :       pari_err_VAR(f, x,relpol);
    1004             :   }
    1005       24570 :   retmkpolmod(scalarpol(x, varn(relpol)), RgX_copy(relpol));
    1006             : }
    1007             : 
    1008             : GEN
    1009        3129 : matbasistoalg(GEN nf,GEN x)
    1010             : {
    1011             :   long i, j, li, lx;
    1012        3129 :   GEN z = cgetg_copy(x, &lx);
    1013             : 
    1014        3129 :   if (lx == 1) return z;
    1015        3115 :   switch(typ(x))
    1016             :   {
    1017             :     case t_VEC: case t_COL:
    1018          28 :       for (i=1; i<lx; i++) gel(z,i) = basistoalg(nf, gel(x,i));
    1019          28 :       return z;
    1020        3087 :     case t_MAT: break;
    1021           0 :     default: pari_err_TYPE("matbasistoalg",x);
    1022             :   }
    1023        3087 :   li = lgcols(x);
    1024       13790 :   for (j=1; j<lx; j++)
    1025             :   {
    1026       10703 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1027       10703 :     gel(z,j) = c;
    1028       10703 :     for (i=1; i<li; i++) gel(c,i) = basistoalg(nf, gel(xj,i));
    1029             :   }
    1030        3087 :   return z;
    1031             : }
    1032             : 
    1033             : GEN
    1034        2106 : matalgtobasis(GEN nf,GEN x)
    1035             : {
    1036             :   long i, j, li, lx;
    1037        2106 :   GEN z = cgetg_copy(x, &lx);
    1038             : 
    1039        2106 :   if (lx == 1) return z;
    1040        2071 :   switch(typ(x))
    1041             :   {
    1042             :     case t_VEC: case t_COL:
    1043        2064 :       for (i=1; i<lx; i++) gel(z,i) = algtobasis(nf, gel(x,i));
    1044        2064 :       return z;
    1045           7 :     case t_MAT: break;
    1046           0 :     default: pari_err_TYPE("matalgtobasis",x);
    1047             :   }
    1048           7 :   li = lgcols(x);
    1049          14 :   for (j=1; j<lx; j++)
    1050             :   {
    1051           7 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1052           7 :     gel(z,j) = c;
    1053           7 :     for (i=1; i<li; i++) gel(c,i) = algtobasis(nf, gel(xj,i));
    1054             :   }
    1055           7 :   return z;
    1056             : }
    1057             : GEN
    1058        3815 : RgM_to_nfM(GEN nf,GEN x)
    1059             : {
    1060             :   long i, j, li, lx;
    1061        3815 :   GEN z = cgetg_copy(x, &lx);
    1062             : 
    1063        3815 :   if (lx == 1) return z;
    1064        3815 :   li = lgcols(x);
    1065       26691 :   for (j=1; j<lx; j++)
    1066             :   {
    1067       22876 :     GEN c = cgetg(li,t_COL), xj = gel(x,j);
    1068       22876 :     gel(z,j) = c;
    1069       22876 :     for (i=1; i<li; i++) gel(c,i) = nf_to_scalar_or_basis(nf, gel(xj,i));
    1070             :   }
    1071        3815 :   return z;
    1072             : }
    1073             : GEN
    1074       38016 : RgC_to_nfC(GEN nf,GEN x)
    1075             : {
    1076       38016 :   long i, lx = lg(x);
    1077       38016 :   GEN z = cgetg(lx, t_COL);
    1078       38016 :   for (i=1; i<lx; i++) gel(z,i) = nf_to_scalar_or_basis(nf, gel(x,i));
    1079       38016 :   return z;
    1080             : }
    1081             : 
    1082             : /* x a t_POLMOD, supposedly in rnf = K[z]/(T), K = Q[y]/(Tnf) */
    1083             : GEN
    1084       60529 : polmod_nffix(const char *f, GEN rnf, GEN x, int lift)
    1085       60529 : { return polmod_nffix2(f, rnf_get_nfpol(rnf), rnf_get_pol(rnf), x,lift); }
    1086             : GEN
    1087       60620 : polmod_nffix2(const char *f, GEN T, GEN relpol, GEN x, int lift)
    1088             : {
    1089       60620 :   if (RgX_equal_var(gel(x,1),relpol))
    1090             :   {
    1091       55370 :     x = gel(x,2);
    1092       55370 :     if (typ(x) == t_POL && varn(x) == varn(relpol))
    1093             :     {
    1094       39655 :       x = RgX_nffix(f, T, x, lift);
    1095       39655 :       switch(lg(x))
    1096             :       {
    1097       11130 :         case 2: return gen_0;
    1098        4193 :         case 3: return gel(x,2);
    1099             :       }
    1100       24332 :       return x;
    1101             :     }
    1102             :   }
    1103       20965 :   return Rg_nffix(f, T, x, lift);
    1104             : }
    1105             : GEN
    1106        1176 : rnfalgtobasis(GEN rnf,GEN x)
    1107             : {
    1108        1176 :   const char *f = "rnfalgtobasis";
    1109        1176 :   pari_sp av = avma;
    1110             :   GEN T, relpol;
    1111             : 
    1112        1176 :   checkrnf(rnf);
    1113        1176 :   relpol = rnf_get_pol(rnf);
    1114        1176 :   T = rnf_get_nfpol(rnf);
    1115        1176 :   switch(typ(x))
    1116             :   {
    1117             :     case t_COL:
    1118          49 :       if (lg(x)-1 != rnf_get_degree(rnf)) pari_err_DIM(f);
    1119          28 :       x = RgV_nffix(f, T, x, 0);
    1120          21 :       return gerepilecopy(av, x);
    1121             : 
    1122             :     case t_POLMOD:
    1123        1043 :       x = polmod_nffix(f, rnf, x, 0);
    1124        1001 :       if (typ(x) != t_POL) break;
    1125         707 :       return gerepileupto(av, mulmat_pol(rnf_get_invzk(rnf), x));
    1126             :     case t_POL:
    1127          56 :       if (varn(x) == varn(T)) { RgX_check_QX(x,f); x = mkpolmod(x,T); break; }
    1128          35 :       x = RgX_nffix(f, T, x, 0);
    1129          28 :       if (degpol(x) >= degpol(relpol)) x = RgX_rem(x,relpol);
    1130          28 :       return gerepileupto(av, mulmat_pol(rnf_get_invzk(rnf), x));
    1131             :   }
    1132         336 :   return gerepileupto(av, scalarcol(x, rnf_get_degree(rnf)));
    1133             : }
    1134             : 
    1135             : /* Given a and b in nf, gives an algebraic integer y in nf such that a-b.y
    1136             :  * is "small" */
    1137             : GEN
    1138         259 : nfdiveuc(GEN nf, GEN a, GEN b)
    1139             : {
    1140         259 :   pari_sp av = avma;
    1141         259 :   a = nfdiv(nf,a,b);
    1142         259 :   return gerepileupto(av, ground(a));
    1143             : }
    1144             : 
    1145             : /* Given a and b in nf, gives a "small" algebraic integer r in nf
    1146             :  * of the form a-b.y */
    1147             : GEN
    1148         259 : nfmod(GEN nf, GEN a, GEN b)
    1149             : {
    1150         259 :   pari_sp av = avma;
    1151         259 :   GEN p1 = gneg_i(nfmul(nf,b,ground(nfdiv(nf,a,b))));
    1152         259 :   return gerepileupto(av, nfadd(nf,a,p1));
    1153             : }
    1154             : 
    1155             : /* Given a and b in nf, gives a two-component vector [y,r] in nf such
    1156             :  * that r=a-b.y is "small". */
    1157             : GEN
    1158         259 : nfdivrem(GEN nf, GEN a, GEN b)
    1159             : {
    1160         259 :   pari_sp av = avma;
    1161         259 :   GEN p1,z, y = ground(nfdiv(nf,a,b));
    1162             : 
    1163         259 :   p1 = gneg_i(nfmul(nf,b,y));
    1164         259 :   z = cgetg(3,t_VEC);
    1165         259 :   gel(z,1) = gcopy(y);
    1166         259 :   gel(z,2) = nfadd(nf,a,p1); return gerepileupto(av, z);
    1167             : }
    1168             : 
    1169             : /*************************************************************************/
    1170             : /**                                                                     **/
    1171             : /**                           (Z_K/I)^*                                 **/
    1172             : /**                                                                     **/
    1173             : /*************************************************************************/
    1174             : /* return sign(sigma_k(x)), x t_COL (integral, primitive) */
    1175             : static long
    1176      372402 : eval_sign(GEN M, GEN x, long k)
    1177             : {
    1178      372402 :   long i, l = lg(x);
    1179      372402 :   GEN z = gel(x,1); /* times M[k,1], which is 1 */
    1180      372402 :   for (i = 2; i < l; i++) z = mpadd(z, mpmul(gcoeff(M,k,i), gel(x,i)));
    1181      372402 :   if (realprec(z) < DEFAULTPREC) pari_err_PREC("nfsign_arch");
    1182      372402 :   return signe(z);
    1183             : }
    1184             : 
    1185             : /* sigma_k(x), assuming x not rational (or nf != Q) */
    1186             : static GEN
    1187         287 : nfembed_i(GEN nf, GEN x, long k)
    1188             : {
    1189             :   long i, l;
    1190             :   GEN z, M;
    1191         287 :   M = nf_get_M(nf); l = lg(M); /* > 2 */
    1192         287 :   z = gel(x,1);
    1193         287 :   for (i=2; i<l; i++) z = gadd(z, gmul(gcoeff(M,k,i), gel(x,i)));
    1194         287 :   return z;
    1195             : }
    1196             : GEN
    1197        1568 : nfembed(GEN nf, GEN x, long k)
    1198             : {
    1199        1568 :   pari_sp av = avma;
    1200        1568 :   nf = checknf(nf);
    1201        1568 :   x = nf_to_scalar_or_basis(nf,x);
    1202        1568 :   if (typ(x) != t_COL) return gerepilecopy(av, x);
    1203           0 :   return gerepileupto(av, nfembed_i(nf,x,k));
    1204             : }
    1205             : 
    1206             : /* pl : requested signs for real embeddings, 0 = no sign constraint */
    1207             : /* FIXME: not rigorous */
    1208             : int
    1209        1113 : nfchecksigns(GEN nf, GEN x, GEN pl)
    1210             : {
    1211        1113 :   pari_sp av = avma;
    1212        1113 :   long l = lg(pl), i;
    1213        1113 :   nf = checknf(nf);
    1214        1113 :   x = nf_to_scalar_or_basis(nf,x);
    1215        1113 :   if (typ(x) != t_COL)
    1216             :   {
    1217         791 :     long s = gsigne(x);
    1218        1631 :     for (i = 1; i < l; i++)
    1219        1036 :       if (pl[i] && pl[i] != s) { avma = av; return 0; }
    1220             :   }
    1221             :   else
    1222             :   {
    1223         525 :     for (i = 1; i < l; i++)
    1224         322 :       if (pl[i] && pl[i] != gsigne(nfembed_i(nf,x,i))) { avma = av; return 0; }
    1225             :   }
    1226         798 :   avma = av; return 1;
    1227             : }
    1228             : 
    1229             : GEN
    1230         595 : vecsmall01_to_indices(GEN v)
    1231             : {
    1232         595 :   long i, k, l = lg(v);
    1233         595 :   GEN p = new_chunk(l) + l;
    1234        1526 :   for (k=1, i=l-1; i; i--)
    1235         931 :     if (v[i]) { *--p = i; k++; }
    1236         595 :   *--p = evallg(k) | evaltyp(t_VECSMALL);
    1237         595 :   avma = (pari_sp)p; return p;
    1238             : }
    1239             : GEN
    1240      734527 : vec01_to_indices(GEN v)
    1241             : {
    1242             :   long i, k, l;
    1243             :   GEN p;
    1244             : 
    1245      734527 :   switch (typ(v))
    1246             :   {
    1247      442118 :    case t_VECSMALL: return v;
    1248      292409 :    case t_VEC: break;
    1249           0 :    default: pari_err_TYPE("vec01_to_indices",v);
    1250             :   }
    1251      292409 :   l = lg(v);
    1252      292409 :   p = new_chunk(l) + l;
    1253      948389 :   for (k=1, i=l-1; i; i--)
    1254      655980 :     if (signe(gel(v,i))) { *--p = i; k++; }
    1255      292409 :   *--p = evallg(k) | evaltyp(t_VECSMALL);
    1256      292409 :   avma = (pari_sp)p; return p;
    1257             : }
    1258             : GEN
    1259        4368 : indices_to_vec01(GEN p, long r)
    1260             : {
    1261        4368 :   long i, l = lg(p);
    1262        4368 :   GEN v = zerovec(r);
    1263        4368 :   for (i = 1; i < l; i++) gel(v, p[i]) = gen_1;
    1264        4368 :   return v;
    1265             : }
    1266             : 
    1267             : /* return (column) vector of R1 signatures of x (0 or 1) */
    1268             : GEN
    1269      525723 : nfsign_arch(GEN nf, GEN x, GEN arch)
    1270             : {
    1271      525723 :   GEN M, V, archp = vec01_to_indices(arch);
    1272      525723 :   long i, s, n = lg(archp)-1;
    1273             :   pari_sp av;
    1274             : 
    1275      525723 :   if (!n) return cgetg(1,t_VECSMALL);
    1276      417468 :   nf = checknf(nf);
    1277      417468 :   if (typ(x) == t_MAT)
    1278             :   { /* factorisation */
    1279      102015 :     GEN g = gel(x,1), e = gel(x,2);
    1280      102015 :     V = zero_zv(n);
    1281      313003 :     for (i=1; i<lg(g); i++)
    1282      210988 :       if (mpodd(gel(e,i)))
    1283      173079 :         Flv_add_inplace(V, nfsign_arch(nf,gel(g,i),archp), 2);
    1284      102015 :     avma = (pari_sp)V; return V;
    1285             :   }
    1286      315453 :   av = avma; V = cgetg(n+1,t_VECSMALL);
    1287      315453 :   x = nf_to_scalar_or_basis(nf, x);
    1288      315453 :   switch(typ(x))
    1289             :   {
    1290             :     case t_INT:
    1291       83408 :       s = signe(x);
    1292       83408 :       if (!s) pari_err_DOMAIN("nfsign_arch","element","=",gen_0,x);
    1293       83408 :       avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
    1294             :     case t_FRAC:
    1295         476 :       s = signe(gel(x,1));
    1296         476 :       avma = av; return const_vecsmall(n, (s < 0)? 1: 0);
    1297             :   }
    1298      231569 :   x = Q_primpart(x); M = nf_get_M(nf);
    1299      231569 :   for (i = 1; i <= n; i++) V[i] = (eval_sign(M, x, archp[i]) < 0)? 1: 0;
    1300      231569 :   avma = (pari_sp)V; return V;
    1301             : }
    1302             : 
    1303             : /* return the vector of signs of x; the matrix of such if x is a vector
    1304             :  * of nf elements */
    1305             : GEN
    1306         378 : nfsign(GEN nf, GEN x)
    1307             : {
    1308             :   long i, l;
    1309             :   GEN arch, S;
    1310             : 
    1311         378 :   nf = checknf(nf);
    1312         378 :   arch = identity_perm( nf_get_r1(nf) );
    1313         378 :   if (typ(x) != t_VEC) return nfsign_arch(nf, x, arch);
    1314          70 :   l = lg(x); S = cgetg(l, t_MAT);
    1315          70 :   for (i=1; i<l; i++) gel(S,i) = nfsign_arch(nf, gel(x,i), arch);
    1316          70 :   return S;
    1317             : }
    1318             : 
    1319             : /* multiply y by t = 1 mod^* f such that sign(x) = sign(y) at arch = divisor[2].
    1320             :  * If x == NULL, make y >> 0 at sarch */
    1321             : GEN
    1322       83384 : set_sign_mod_divisor(GEN nf, GEN x, GEN y, GEN divisor, GEN sarch)
    1323             : {
    1324             :   GEN s, archp, gen;
    1325             :   long nba,i;
    1326       83384 :   if (!sarch) return y;
    1327       83384 :   gen = gel(sarch,2); nba = lg(gen);
    1328       83384 :   if (nba == 1) return y;
    1329             : 
    1330       71659 :   archp = vec01_to_indices(gel(divisor,2));
    1331       71659 :   y = nf_to_scalar_or_basis(nf, y);
    1332       71659 :   s = nfsign_arch(nf, y, archp);
    1333       71659 :   if (x) Flv_add_inplace(s, nfsign_arch(nf, x, archp), 2);
    1334       71659 :   s = Flm_Flc_mul(gel(sarch,3), s, 2);
    1335      168651 :   for (i=1; i<nba; i++)
    1336       96992 :     if (s[i]) y = nfmul(nf,y,gel(gen,i));
    1337       71659 :   if (typ(y) != t_COL) y = scalarcol_shallow(y, nf_get_degree(nf));
    1338       71659 :   return y;
    1339             : }
    1340             : 
    1341             : /* given an element x in Z_K and an integral ideal y in HNF, coprime with x,
    1342             :    outputs an element inverse of x modulo y */
    1343             : GEN
    1344         378 : nfinvmodideal(GEN nf, GEN x, GEN y)
    1345             : {
    1346         378 :   pari_sp av = avma;
    1347         378 :   GEN a, yZ = gcoeff(y,1,1);
    1348             : 
    1349         378 :   if (is_pm1(yZ)) return zerocol( nf_get_degree(nf) );
    1350         378 :   x = nf_to_scalar_or_basis(nf, x);
    1351         378 :   if (typ(x) == t_INT) return gerepileupto(av, Fp_inv(x, yZ));
    1352             : 
    1353         231 :   a = hnfmerge_get_1(idealhnf_principal(nf,x), y);
    1354         231 :   if (!a) pari_err_INV("nfinvmodideal", x);
    1355         231 :   return gerepileupto(av, ZC_hnfrem(nfdiv(nf,a,x), y));
    1356             : }
    1357             : 
    1358             : static GEN
    1359      282737 : nfsqrmodideal(GEN nf, GEN x, GEN id) {
    1360      282737 :   return ZC_hnfrem(nfsqri(nf,x), id);
    1361             : }
    1362             : static GEN
    1363      601147 : nfmulmodideal(GEN nf, GEN x, GEN y, GEN id) {
    1364      601147 :   return x? ZC_hnfrem(nfmuli(nf,x,y), id): y;
    1365             : }
    1366             : /* assume x integral, k integer, A in HNF */
    1367             : GEN
    1368      384460 : nfpowmodideal(GEN nf,GEN x,GEN k,GEN A)
    1369             : {
    1370      384460 :   long s = signe(k);
    1371             :   pari_sp av;
    1372             :   GEN y;
    1373             : 
    1374      384460 :   if (!s) return gen_1;
    1375      384460 :   av = avma;
    1376      384460 :   x = nf_to_scalar_or_basis(nf, x);
    1377      384460 :   if (typ(x) != t_COL) return Fp_pow(x, k, gcoeff(A,1,1));
    1378      199549 :   if (s < 0) { x = nfinvmodideal(nf, x,A); k = absi(k); }
    1379      199549 :   for(y = NULL;;)
    1380             :   {
    1381      482286 :     if (mpodd(k)) y = nfmulmodideal(nf,y,x,A);
    1382      482286 :     k = shifti(k,-1); if (!signe(k)) break;
    1383      282737 :     x = nfsqrmodideal(nf,x,A);
    1384      282737 :   }
    1385      199549 :   return gerepileupto(av, y);
    1386             : }
    1387             : 
    1388             : /* a * g^n mod id */
    1389             : static GEN
    1390      325478 : elt_mulpow_modideal(GEN nf, GEN a, GEN g, GEN n, GEN id)
    1391             : {
    1392      325478 :   return nfmulmodideal(nf, a, nfpowmodideal(nf,g,n,id), id);
    1393             : }
    1394             : 
    1395             : /* assume (num(g[i]), id) = 1 for all i. Return prod g[i]^e[i] mod id.
    1396             :  * EX = multiple of exponent of (O_K/id)^* */
    1397             : GEN
    1398      112315 : famat_to_nf_modideal_coprime(GEN nf, GEN g, GEN e, GEN id, GEN EX)
    1399             : {
    1400      112315 :   GEN plus = NULL, minus = NULL, idZ = gcoeff(id,1,1);
    1401      112315 :   pari_sp av = avma;
    1402      112315 :   long i, lx = lg(g);
    1403      112315 :   GEN EXo2 = (expi(EX) > 10)? shifti(EX,-1): NULL;
    1404             : 
    1405      112315 :   if (is_pm1(idZ)) lx = 1; /* id = Z_K */
    1406      412419 :   for (i=1; i<lx; i++)
    1407             :   {
    1408      300104 :     GEN h, n = centermodii(gel(e,i), EX, EXo2);
    1409      300104 :     long sn = signe(n);
    1410      300104 :     if (!sn) continue;
    1411             : 
    1412      245969 :     h = nf_to_scalar_or_basis(nf, gel(g,i));
    1413      245969 :     switch(typ(h))
    1414             :     {
    1415      160191 :       case t_INT: break;
    1416             :       case t_FRAC:
    1417           0 :         h = Fp_div(gel(h,1), gel(h,2), idZ); break;
    1418             :       default:
    1419             :       {
    1420             :         GEN dh;
    1421       85778 :         h = Q_remove_denom(h, &dh);
    1422       85778 :         if (dh) h = FpC_Fp_mul(h, Fp_inv(dh,idZ), idZ);
    1423             :       }
    1424             :     }
    1425      245969 :     if (sn > 0)
    1426      244310 :       plus = elt_mulpow_modideal(nf, plus, h, n, id);
    1427             :     else /* sn < 0 */
    1428        1659 :       minus = elt_mulpow_modideal(nf, minus, h, absi(n), id);
    1429             : 
    1430      245969 :     if (gc_needed(av, 2))
    1431             :     {
    1432           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"famat_to_nf_modideal_coprime");
    1433           0 :       if (!plus) plus = gen_0;
    1434           0 :       if (!minus) minus = gen_0;
    1435           0 :       gerepileall(av,2, &plus, &minus);
    1436           0 :       if (isintzero(plus)) plus = NULL;
    1437           0 :       if (isintzero(minus)) minus = NULL;
    1438             :     }
    1439             :   }
    1440      112315 :   if (minus)
    1441         378 :     plus = nfmulmodideal(nf, plus, nfinvmodideal(nf,minus,id), id);
    1442      112315 :   return plus? plus: scalarcol_shallow(gen_1, lg(id)-1);
    1443             : }
    1444             : 
    1445             : /* given 2 integral ideals x, y in HNF s.t x | y | x^2, compute the quotient
    1446             :    (1+x)/(1+y) in the form [[cyc],[gen]], if U != NULL, set *U := ux^-1 */
    1447             : static GEN
    1448       12985 : zidealij(GEN x, GEN y, GEN *U)
    1449             : {
    1450             :   GEN G, cyc;
    1451             :   long j, N;
    1452             : 
    1453             :   /* x^(-1) y = relations between the 1 + x_i (HNF) */
    1454       12985 :   cyc = ZM_snf_group(hnf_solve(x, y), U, &G);
    1455       12985 :   N = lg(cyc); G = ZM_mul(x,G); settyp(G, t_VEC); /* new generators */
    1456       60438 :   for (j=1; j<N; j++)
    1457             :   {
    1458       47453 :     GEN c = gel(G,j);
    1459       47453 :     gel(c,1) = addiu(gel(c,1), 1); /* 1 + g_j */
    1460       47453 :     if (ZV_isscalar(c)) gel(G,j) = gel(c,1);
    1461             :   }
    1462       12985 :   if (U) *U = RgM_mul(*U, RgM_inv(x));
    1463       12985 :   return mkvec2(cyc, G);
    1464             : }
    1465             : 
    1466             : static GEN
    1467      299031 : Fq_FpXQ_log(GEN a, GEN g, GEN ord, GEN T, GEN p)
    1468             : {
    1469      299031 :   if (!T) return Fp_log(a,g,ord,p);
    1470      120751 :   if (typ(a)==t_INT) return Fp_FpXQ_log(a,g,ord,T,p);
    1471       83777 :   return FpXQ_log(a,g,ord,T,p);
    1472             : }
    1473             : /* same in nf.zk / pr */
    1474             : static GEN
    1475      299031 : nf_log(GEN nf, GEN a, GEN g, GEN ord, GEN pr)
    1476             : {
    1477      299031 :   pari_sp av = avma;
    1478      299031 :   GEN T,p, modpr = nf_to_Fq_init(nf, &pr, &T, &p);
    1479      299031 :   GEN A = nf_to_Fq(nf,a,modpr);
    1480      299031 :   GEN G = nf_to_Fq(nf,g,modpr);
    1481      299031 :   return gerepileuptoint(av, Fq_FpXQ_log(A,G,ord,T,p));
    1482             : }
    1483             : 
    1484             : /* Product of cyc entries, with cyc = diagonal(Smith Normal Form), assumed != 0.
    1485             :  * Set L to the index of the last non-trivial (!= 1) entry */
    1486             : GEN
    1487        7434 : detcyc(GEN cyc, long *L)
    1488             : {
    1489        7434 :   pari_sp av = avma;
    1490        7434 :   long i, l = lg(cyc);
    1491        7434 :   GEN s = gel(cyc,1);
    1492             : 
    1493        7434 :   if (l == 1) { *L = 1; return gen_1; }
    1494        9751 :   for (i=2; i<l; i++)
    1495             :   {
    1496        2940 :     GEN t = gel(cyc,i);
    1497        2940 :     if (is_pm1(t)) break;
    1498        2359 :     s = mulii(s, t);
    1499             :   }
    1500        7392 :   *L = i; return i <= 2? icopy(s): gerepileuptoint(av,s);
    1501             : }
    1502             : 
    1503             : 
    1504             : /* lg(x) > 1, x + 1; shallow */
    1505             : static GEN
    1506        3500 : ZC_add1(GEN x)
    1507             : {
    1508        3500 :   long i, l = lg(x);
    1509        3500 :   GEN y = cgetg(l, t_COL);
    1510        3500 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    1511        3500 :   gel(y,1) = addiu(gel(x,1), 1); return y;
    1512             : }
    1513             : /* lg(x) > 1, x - 1; shallow */
    1514             : static GEN
    1515        1561 : ZC_sub1(GEN x)
    1516             : {
    1517        1561 :   long i, l = lg(x);
    1518        1561 :   GEN y = cgetg(l, t_COL);
    1519        1561 :   for (i = 2; i < l; i++) gel(y,i) = gel(x,i);
    1520        1561 :   gel(y,1) = subiu(gel(x,1), 1); return y;
    1521             : }
    1522             : 
    1523             : /* x,y are t_INT or ZC */
    1524             : static GEN
    1525           0 : zkadd(GEN x, GEN y)
    1526             : {
    1527           0 :   long tx = typ(x);
    1528           0 :   if (tx == typ(y))
    1529           0 :     return tx == t_INT? addii(x,y): ZC_add(x,y);
    1530             :   else
    1531           0 :     return tx == t_INT? ZC_Z_add(y,x): ZC_Z_add(x,y);
    1532             : }
    1533             : /* x a t_INT or ZC, x+1; shallow */
    1534             : static GEN
    1535        4235 : zkadd1(GEN x)
    1536             : {
    1537        4235 :   long tx = typ(x);
    1538        4235 :   return tx == t_INT? addiu(x,1): ZC_add1(x);
    1539             : }
    1540             : /* x a t_INT or ZC, x-1; shallow */
    1541             : static GEN
    1542        4235 : zksub1(GEN x)
    1543             : {
    1544        4235 :   long tx = typ(x);
    1545        4235 :   return tx == t_INT? subiu(x,1): ZC_sub1(x);
    1546             : }
    1547             : /* x,y are t_INT or ZC; x - y */
    1548             : static GEN
    1549           0 : zksub(GEN x, GEN y)
    1550             : {
    1551           0 :   long tx = typ(x), ty = typ(y);
    1552           0 :   if (tx == ty)
    1553           0 :     return tx == t_INT? subii(x,y): ZC_sub(x,y);
    1554             :   else
    1555           0 :     return tx == t_INT? Z_ZC_sub(x,y): ZC_Z_sub(x,y);
    1556             : }
    1557             : /* x is t_INT or ZM (mult. map), y is t_INT or ZC; x * y */
    1558             : static GEN
    1559        4235 : zkmul(GEN x, GEN y)
    1560             : {
    1561        4235 :   long tx = typ(x), ty = typ(y);
    1562        4235 :   if (ty == t_INT)
    1563        2674 :     return tx == t_INT? mulii(x,y): ZC_Z_mul(gel(x,1),y);
    1564             :   else
    1565        1561 :     return tx == t_INT? ZC_Z_mul(y,x): ZM_ZC_mul(x,y);
    1566             : }
    1567             : 
    1568             : /* (U,V) = 1 coprime ideals. Want z = x mod U, = y mod V; namely
    1569             :  * z =vx + uy = v(x-y) + y, where u + v = 1, u in U, v in V.
    1570             :  * zkc = [v, UV], v a t_INT or ZM (mult. by v map), UV a ZM (ideal in HNF);
    1571             :  * shallow */
    1572             : GEN
    1573           0 : zkchinese(GEN zkc, GEN x, GEN y)
    1574             : {
    1575           0 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd(zkmul(v, zksub(x,y)), y);
    1576           0 :   return (typ(z) == t_INT)? modii(z, gcoeff(UV,1,1)): ZC_hnfrem(z, UV);
    1577             : }
    1578             : /* special case z = x mod U, = 1 mod V; shallow */
    1579             : GEN
    1580        4235 : zkchinese1(GEN zkc, GEN x)
    1581             : {
    1582        4235 :   GEN v = gel(zkc,1), UV = gel(zkc,2), z = zkadd1(zkmul(v, zksub1(x)));
    1583        4235 :   return (typ(z) == t_INT)? z: ZC_hnfrem(z, UV);
    1584             : }
    1585             : static GEN
    1586         518 : zkVchinese1(GEN zkc, GEN v)
    1587             : {
    1588             :   long i, ly;
    1589         518 :   GEN y = cgetg_copy(v, &ly);
    1590         518 :   for (i=1; i<ly; i++) gel(y,i) = zkchinese1(zkc, gel(v,i));
    1591         518 :   return y;
    1592             : }
    1593             : 
    1594             : /* prepare to solve z = x (mod A), z = y mod (B) [zkchinese or zkchinese1] */
    1595             : GEN
    1596        2478 : zkchineseinit(GEN nf, GEN A, GEN B, GEN AB)
    1597             : {
    1598             :   GEN v;
    1599        2478 :   nf = checknf(nf);
    1600        2478 :   v = idealaddtoone_i(nf, A, B);
    1601        2471 :   return mkvec2(zk_scalar_or_multable(nf,v), AB);
    1602             : }
    1603             : /* prepare to solve z = x (mod A), z = 1 mod (B)
    1604             :  * and then         z = 1 (mod A), z = y mod (B) [zkchinese1 twice] */
    1605             : static GEN
    1606         259 : zkchinese1init2(GEN nf, GEN A, GEN B, GEN AB)
    1607             : {
    1608         259 :   GEN zkc = zkchineseinit(nf, A, B, AB);
    1609         259 :   GEN mv = gel(zkc,1), mu;
    1610         259 :   if (typ(mv) == t_INT) return mkvec2(zkc, mkvec2(subui(1,mv),AB));
    1611         238 :   mu = RgM_Rg_add_shallow(ZM_neg(mv), gen_1);
    1612         238 :   return mkvec2(mkvec2(mv,AB), mkvec2(mu,AB));
    1613             : }
    1614             : 
    1615             : /* Given an ideal pr^ep, and an integral ideal x (in HNF form) compute a list
    1616             :  * of vectors,corresponding to the abelian groups (O_K/pr)^*, and
    1617             :  * 1 + pr^i/ 1 + pr^min(2i, ep), i = 1,...
    1618             :  * Each vector has 5 components as follows :
    1619             :  * [[cyc],[g],[g'],[sign],U.X^-1].
    1620             :  * cyc   = type as abelian group
    1621             :  * g, g' = generators. (g',x) = 1, not necessarily so for g
    1622             :  * sign  = vector of the sign(g') at arch.
    1623             :  * If x = NULL, the original ideal was a prime power */
    1624             : static GEN
    1625       12922 : zprimestar(GEN nf, GEN pr, GEN ep, GEN x, GEN arch)
    1626             : {
    1627       12922 :   pari_sp av = avma;
    1628       12922 :   long a, e = itos(ep), f = pr_get_f(pr);
    1629       12922 :   GEN p = pr_get_p(pr), list, g, g0, y, uv, prb, pre;
    1630             :   ulong mask;
    1631             : 
    1632       12922 :   if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",e,pr);
    1633       12922 :   if (f == 1)
    1634        6279 :     g = pgener_Fp(p);
    1635             :   else
    1636             :   {
    1637        6643 :     GEN T, modpr = zk_to_Fq_init(nf, &pr, &T, &p);
    1638        6643 :     g = Fq_to_nf(gener_FpXQ(T,p,NULL), modpr);
    1639        6643 :     g = poltobasis(nf, g);
    1640             :   }
    1641             :   /* g generates  (Z_K / pr)^* */
    1642       12922 :   prb = idealhnf_two(nf,pr);
    1643       12922 :   pre = (e==1)? prb: idealpow(nf,pr,ep);
    1644       12922 :   if (x)
    1645             :   {
    1646        2219 :     uv = zkchineseinit(nf, idealdivpowprime(nf,x,pr,ep), pre, x);
    1647        2212 :     g0 = zkchinese1(uv, g);
    1648             :   }
    1649             :   else
    1650             :   {
    1651       10703 :     uv = NULL; /* gcc -Wall */
    1652       10703 :     g0 = g;
    1653             :   }
    1654             : 
    1655       12915 :   y = mkvec5(mkvec(subiu(powiu(p,f), 1)),
    1656             :              mkvec(g),
    1657             :              mkvec(g0),
    1658             :              mkvec(nfsign_arch(nf,g0,arch)),
    1659             :              gen_1);
    1660       12915 :   if (e == 1) return gerepilecopy(av, mkvec(y));
    1661        7784 :   list = vectrunc_init(e+1);
    1662        7784 :   vectrunc_append(list, y);
    1663        7784 :   mask = quadratic_prec_mask(e);
    1664        7784 :   a = 1;
    1665       27335 :   while (mask > 1)
    1666             :   {
    1667       11767 :     GEN pra = prb, gen, z, s, U;
    1668       11767 :     long i, l, b = a << 1;
    1669             : 
    1670       11767 :     if (mask & 1) b--;
    1671       11767 :     mask >>= 1;
    1672             :     /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
    1673       11767 :     if(DEBUGLEVEL>3) err_printf("  treating a = %ld, b = %ld\n",a,b);
    1674       11767 :     prb = (b >= e)? pre: idealpows(nf,pr,b);
    1675       11767 :     z = zidealij(pra, prb, &U);
    1676       11767 :     gen = leafcopy(gel(z,2));
    1677       11767 :     s = cgetg_copy(gen, &l);
    1678       57778 :     for (i = 1; i < l; i++)
    1679             :     {
    1680       46011 :       if (x) gel(gen,i) = zkchinese1(uv, gel(gen,i));
    1681       46011 :       gel(s,i) = nfsign_arch(nf, gel(gen,i), arch);
    1682             :     }
    1683       11767 :     y = mkvec5(gel(z,1), gel(z,2), gen, s, U);
    1684       11767 :     vectrunc_append(list, y);
    1685       11767 :     a = b;
    1686             :   }
    1687        7784 :   return gerepilecopy(av, list);
    1688             : }
    1689             : 
    1690             : /* increment y, which runs through [-d,d]^k. Return 0 when done. */
    1691             : static int
    1692        3192 : increment(GEN y, long k, long d)
    1693             : {
    1694        3192 :   long i = 0, j;
    1695             :   do
    1696             :   {
    1697        4480 :     if (++i > k) return 0;
    1698        4480 :     y[i]++;
    1699        4480 :   } while (y[i] > d);
    1700        3192 :   for (j = 1; j < i; j++) y[j] = -d;
    1701        3192 :   return 1;
    1702             : }
    1703             : 
    1704             : GEN
    1705        1729 : archstar_full_rk(GEN x, GEN bas, GEN v, GEN gen)
    1706             : {
    1707        1729 :   long i, r, lgmat, N = lg(bas)-1, nba = nbrows(v);
    1708        1729 :   GEN lambda = cgetg(N+1, t_VECSMALL), mat = cgetg(nba+1,t_MAT);
    1709             : 
    1710        1729 :   lgmat = lg(v); setlg(mat, lgmat+1);
    1711        1729 :   for (i = 1; i < lgmat; i++) gel(mat,i) = gel(v,i);
    1712        1729 :   for (     ; i <= nba; i++)  gel(mat,i) = cgetg(nba+1, t_VECSMALL);
    1713             : 
    1714        1729 :   if (x) { x = ZM_lll(x, 0.75, LLL_INPLACE); bas = RgV_RgM_mul(bas, x); }
    1715             : 
    1716        1729 :   for (r=1;; r++)
    1717             :   { /* reset */
    1718        1729 :     (void)vec_setconst(lambda, (GEN)0);
    1719        1729 :     if (!x) lambda[1] = r;
    1720        4921 :     while (increment(lambda, N, r))
    1721             :     {
    1722        3192 :       pari_sp av1 = avma;
    1723        3192 :       GEN a = RgM_zc_mul(bas, lambda), c = gel(mat,lgmat);
    1724       10122 :       for (i = 1; i <= nba; i++)
    1725             :       {
    1726        6930 :         GEN t = x? gadd(gel(a,i), gen_1): gel(a,i);
    1727        6930 :         c[i] = (gsigne(t) < 0)? 1: 0;
    1728             :       }
    1729        3192 :       avma = av1; if (Flm_deplin(mat, 2)) continue;
    1730             : 
    1731             :       /* c independent of previous sign vectors */
    1732        1764 :       if (!x) a = zc_to_ZC(lambda);
    1733             :       else
    1734             :       {
    1735        1316 :         a = ZM_zc_mul(x, lambda);
    1736        1316 :         gel(a,1) = addis(gel(a,1), 1);
    1737             :       }
    1738        1764 :       gel(gen,lgmat) = a;
    1739        1764 :       if (lgmat++ == nba) {
    1740        1729 :         mat = Flm_inv(mat,2); /* full rank */
    1741        3458 :         settyp(mat, t_VEC); return mat;
    1742             :       }
    1743          35 :       setlg(mat,lgmat+1);
    1744             :     }
    1745           0 :   }
    1746             : }
    1747             : 
    1748             : /* x non-0 integral ideal in HNF, compute elements in 1+x (in x, if x = zk)
    1749             :  * whose sign matrix is invertible. archp in 'indices' format */
    1750             : GEN
    1751       14091 : nfarchstar(GEN nf, GEN x, GEN archp)
    1752             : {
    1753       14091 :   long nba = lg(archp) - 1;
    1754             :   GEN cyc, gen, mat;
    1755             : 
    1756       14091 :   if (!nba)
    1757        9037 :     cyc = gen = mat = cgetg(1, t_VEC);
    1758             :   else
    1759             :   {
    1760        5054 :     GEN xZ = gcoeff(x,1,1), gZ;
    1761        5054 :     pari_sp av = avma;
    1762        5054 :     if (is_pm1(xZ)) x = NULL; /* x = O_K */
    1763        5054 :     gZ = x? subsi(1, xZ): gen_m1; /* gZ << 0, gZ = 1 mod x */
    1764        5054 :     if (nba == 1)
    1765             :     {
    1766        3332 :       gen = mkvec(gZ);
    1767        3332 :       mat = mkvec( mkvecsmall(1) );
    1768             :     }
    1769             :     else
    1770             :     {
    1771        1722 :       GEN bas = nf_get_M(nf);
    1772        1722 :       if (lgcols(bas) > lg(archp)) bas = rowpermute(bas, archp);
    1773        1722 :       gen = cgetg(nba+1,t_VEC);
    1774        1722 :       gel(gen,1) = gZ;
    1775        1722 :       mat = archstar_full_rk(x, bas, mkmat(const_vecsmall(nba,1)), gen);
    1776        1722 :       gerepileall(av,2,&gen,&mat);
    1777             :     }
    1778        5054 :     cyc = const_vec(nba, gen_2);
    1779             :   }
    1780       14091 :   return mkvec3(cyc,gen,mat);
    1781             : }
    1782             : 
    1783             : static GEN
    1784      220857 : apply_U(GEN U, GEN a)
    1785             : {
    1786             :   GEN e;
    1787      220857 :   if (typ(a) == t_INT)
    1788       79065 :     e = RgC_Rg_mul(gel(U,1), subis(a, 1));
    1789             :   else
    1790             :   { /* t_COL */
    1791      141792 :     GEN t = gel(a,1);
    1792      141792 :     gel(a,1) = addsi(-1, gel(a,1)); /* a -= 1 */
    1793      141792 :     e = RgM_RgC_mul(U, a);
    1794      141792 :     gel(a,1) = t; /* restore */
    1795             :   }
    1796      220857 :   return e;
    1797             : }
    1798             : /* a in Z_K (t_COL or t_INT), pr prime ideal, prk = pr^k,
    1799             :  * list = zprimestar(nf, pr, k, ...)  */
    1800             : static GEN
    1801      299031 : zlog_pk(GEN nf, GEN a, GEN y, GEN pr, GEN prk, GEN list, GEN *psigne)
    1802             : {
    1803      299031 :   long i,j, llist = lg(list)-1;
    1804      818744 :   for (j = 1; j <= llist; j++)
    1805             :   {
    1806      519720 :     GEN L = gel(list,j), e;
    1807      519720 :     GEN cyc = gel(L,1), gen = gel(L,2), s = gel(L,4), U = gel(L,5);
    1808      519720 :     if (j == 1)
    1809      299031 :       e = mkcol( nf_log(nf, a, gel(gen,1), gel(cyc,1), pr) );
    1810             :     else
    1811      220689 :       e = apply_U(U, a);
    1812             :     /* here lg(e) == lg(cyc) */
    1813     1588135 :     for (i = 1; i < lg(cyc); i++)
    1814             :     {
    1815             :       GEN t;
    1816     1068422 :       if (typ(gel(e,i)) != t_INT) pari_err_COPRIME("zlog_pk", a, pr);
    1817     1068415 :       t = modii(negi(gel(e,i)), gel(cyc,i));
    1818     1068415 :       gel(++y,0) = negi(t); if (!signe(t)) continue;
    1819             : 
    1820      336177 :       if (mod2(t)) Flv_add_inplace(*psigne, gel(s,i), 2);
    1821      336177 :       if (j != llist) a = elt_mulpow_modideal(nf, a, gel(gen,i), t, prk);
    1822             :     }
    1823             :   }
    1824      299024 :   return y;
    1825             : }
    1826             : 
    1827             : static void
    1828      287670 : zlog_add_sign(GEN y0, GEN sgn, GEN lists)
    1829             : {
    1830             :   GEN y, s;
    1831             :   long i;
    1832      575340 :   if (!sgn) return;
    1833      287670 :   y = y0 + lg(y0);
    1834      287670 :   s = Flm_Flc_mul(gmael(lists, lg(lists)-1, 3), sgn, 2);
    1835      287670 :   for (i = lg(s)-1; i > 0; i--) gel(--y,0) = s[i]? gen_1: gen_0;
    1836             : }
    1837             : 
    1838             : static GEN
    1839       95344 : famat_zlog(GEN nf, GEN fa, GEN sgn, GEN bid)
    1840             : {
    1841       95344 :   GEN g = gel(fa,1), e = gel(fa,2);
    1842       95344 :   GEN vp = gmael(bid, 3,1), ep = gmael(bid, 3,2);
    1843       95344 :   GEN arch = bid_get_arch(bid);
    1844       95344 :   GEN cyc = bid_get_cyc(bid), lists = gel(bid,4), U = bid_get_U(bid);
    1845       95344 :   GEN y0, x, y, EX = gel(cyc,1);
    1846             :   long i, l;
    1847             : 
    1848       95344 :   y0 = y = cgetg(lg(U), t_COL);
    1849       95344 :   if (!sgn) sgn = nfsign_arch(nf, mkmat2(g,e), arch);
    1850       95344 :   l = lg(vp);
    1851      189480 :   for (i=1; i < l; i++)
    1852             :   {
    1853       94136 :     GEN pr = gel(vp,i), prk, ex;
    1854       94136 :     if (l == 2) {
    1855       66970 :       prk = bid_get_ideal(bid);
    1856       66970 :       ex = EX;
    1857             :     } else { /* try to improve EX: should be group exponent mod prf, not f */
    1858       27166 :       GEN k = gel(ep,i);
    1859       27166 :       prk = idealpow(nf, pr, k);
    1860             :       /* upper bound: gcd(EX, (Nv-1)p^(k-1)) = (Nv-1) p^min(k-1,v_p(EX)) */
    1861       27166 :       ex = subis(pr_norm(pr),1);
    1862       27166 :       if (!is_pm1(k)) {
    1863        3563 :         GEN p = pr_get_p(pr), k_1 = subis(k,1);
    1864        3563 :         long v = Z_pval(EX, p);
    1865        3563 :         if (cmpui(v, k_1) > 0) v = itos(k_1);
    1866        3563 :         if (v) ex = mulii(ex, powiu(p, v));
    1867             :       }
    1868             :     }
    1869       94136 :     x = famat_makecoprime(nf, g, e, pr, prk, ex);
    1870       94136 :     y = zlog_pk(nf, x, y, pr, prk, gel(lists,i), &sgn);
    1871             :   }
    1872       95344 :   zlog_add_sign(y0, sgn, lists);
    1873       95344 :   return y0;
    1874             : }
    1875             : 
    1876             : static GEN
    1877      130208 : get_index(GEN lists)
    1878             : {
    1879      130208 :   long t = 0, j, k, l = lg(lists)-1;
    1880      130208 :   GEN L, ind = cgetg(l+1, t_VECSMALL);
    1881             : 
    1882      287041 :   for (k = 1; k < l; k++)
    1883             :   {
    1884      156833 :     L = gel(lists,k);
    1885      156833 :     ind[k] = t;
    1886      156833 :     for (j=1; j<lg(L); j++) t += lg(gmael(L,j,1)) - 1;
    1887             :   }
    1888             :   /* for arch. components */
    1889      130208 :   ind[k] = t; return ind;
    1890             : }
    1891             : 
    1892             : static void
    1893      130208 : init_zlog(zlog_S *S, long n, GEN P, GEN e, GEN arch, GEN lists, GEN U)
    1894             : {
    1895      130208 :   S->n = n;
    1896      130208 :   S->U = U;
    1897      130208 :   S->P = P;
    1898      130208 :   S->e = e;
    1899      130208 :   S->archp = vec01_to_indices(arch);
    1900      130208 :   S->lists = lists;
    1901      130208 :   S->ind = get_index(lists);
    1902      130208 : }
    1903             : void
    1904      118532 : init_zlog_bid(zlog_S *S, GEN bid)
    1905             : {
    1906      118532 :   GEN fa = bid_get_fact(bid), lists = gel(bid,4), U = bid_get_U(bid);
    1907      118532 :   GEN arch = gel(bid_get_mod(bid), 2);
    1908      118532 :   init_zlog(S, lg(U)-1, gel(fa,1), gel(fa,2), arch, lists, U);
    1909      118532 : }
    1910             : 
    1911             : /* Return decomposition of a on the S->n successive generators contained in
    1912             :  * S->lists. If index !=0, do the computation for the corresponding prime
    1913             :  * ideal and set to 0 the other components. */
    1914             : static GEN
    1915      175596 : zlog_ind(GEN nf, GEN a, zlog_S *S, GEN sgn, long index)
    1916             : {
    1917      175596 :   GEN y0 = zerocol(S->n), y;
    1918      175596 :   pari_sp av = avma;
    1919             :   long k, kmin, kmax;
    1920             : 
    1921      175596 :   a = nf_to_scalar_or_basis(nf,a);
    1922      175596 :   if (index)
    1923             :   {
    1924       58926 :     kmin = kmax = index;
    1925       58926 :     y = y0 + S->ind[index];
    1926             :   }
    1927             :   else
    1928             :   {
    1929      116670 :     kmin = 1; kmax = lg(S->P)-1;
    1930      116670 :     y = y0;
    1931             :   }
    1932      175596 :   if (!sgn) sgn = nfsign_arch(nf, a, S->archp);
    1933      376634 :   for (k = kmin; k <= kmax; k++)
    1934             :   {
    1935      201045 :     GEN list= gel(S->lists,k);
    1936      201045 :     GEN pr  = gel(S->P,k);
    1937      201045 :     GEN prk = idealpow(nf, pr, gel(S->e,k));
    1938      201045 :     y = zlog_pk(nf, a, y, pr, prk, list, &sgn);
    1939             :   }
    1940      175589 :   zlog_add_sign(y0, sgn, S->lists);
    1941      175589 :   return gerepilecopy(av, y0);
    1942             : }
    1943             : /* sgn = sign(a, S->arch) or NULL if unknown */
    1944             : GEN
    1945      116670 : zlog(GEN nf, GEN a, GEN sgn, zlog_S *S) { return zlog_ind(nf, a, S, sgn, 0); }
    1946             : 
    1947             : /* Log on bid.gen of generators of P_{1,I pr^{e-1}} / P_{1,I pr^e} (I,pr) = 1,
    1948             :  * defined implicitly via CRT. 'index' is the index of pr in modulus
    1949             :  * factorization */
    1950             : GEN
    1951        9835 : log_gen_pr(zlog_S *S, long index, GEN nf, long e)
    1952             : {
    1953        9835 :   long i, l, yind = S->ind[index];
    1954        9835 :   GEN y, A, L, L2 = gel(S->lists,index);
    1955             : 
    1956        9835 :   if (e == 1)
    1957             :   {
    1958        6671 :     L = gel(L2,1);
    1959        6671 :     y = col_ei(S->n, yind+1);
    1960        6671 :     zlog_add_sign(y, gmael(L,4,1), S->lists);
    1961        6671 :     retmkmat( ZM_ZC_mul(S->U, y) );
    1962             :   }
    1963             :   else
    1964             :   {
    1965        3164 :     GEN prk, g, pr = gel(S->P,index);
    1966        3164 :     long narchp = lg(S->archp)-1;
    1967             : 
    1968        3164 :     if (e == 2)
    1969        1974 :       L = gel(L2,2);
    1970             :     else
    1971        1190 :       L = zidealij(idealpows(nf,pr,e-1), idealpows(nf,pr,e), NULL);
    1972        3164 :     g = gel(L,2);
    1973        3164 :     l = lg(g); A = cgetg(l, t_MAT);
    1974        3164 :     prk = idealpow(nf, pr, gel(S->e,index));
    1975        7014 :     for (i = 1; i < l; i++)
    1976             :     {
    1977        3850 :       GEN G = gel(g,i), sgn = zero_zv(narchp); /*positive at f_oo*/
    1978        3850 :       y = zerocol(S->n);
    1979        3850 :       (void)zlog_pk(nf, G, y + yind, pr, prk, L2, &sgn);
    1980        3850 :       zlog_add_sign(y, sgn, S->lists);
    1981        3850 :       gel(A,i) = y;
    1982             :     }
    1983        3164 :     return ZM_mul(S->U, A);
    1984             :   }
    1985             : }
    1986             : /* Log on bid.gen of generator of P_{1,f} / P_{1,f v[index]}
    1987             :  * v = vector of r1 real places */
    1988             : GEN
    1989        6216 : log_gen_arch(zlog_S *S, long index)
    1990             : {
    1991        6216 :   GEN y = zerocol(S->n);
    1992        6216 :   zlog_add_sign(y, vecsmall_ei(lg(S->archp)-1, index), S->lists);
    1993        6216 :   return ZM_ZC_mul(S->U, y);
    1994             : }
    1995             : 
    1996             : /* add [h,cyc] or [h,cyc,gen] to bid */
    1997             : static void
    1998       13566 : add_grp(GEN nf, GEN u1, GEN cyc, GEN gen, GEN bid)
    1999             : {
    2000       13566 :   GEN h = ZV_prod(cyc);
    2001       13566 :   if (u1)
    2002             :   {
    2003        7406 :     GEN G = mkvec3(h,cyc,NULL/*dummy, bid[2] needed below*/);
    2004        7406 :     gel(bid,2) = G;
    2005        7406 :     if (u1 != gen_1)
    2006             :     {
    2007        5887 :       long i, c = lg(u1);
    2008        5887 :       GEN g = cgetg(c,t_VEC);
    2009       18249 :       for (i=1; i<c; i++)
    2010       12362 :         gel(g,i) = famat_to_nf_moddivisor(nf, gen, gel(u1,i), bid);
    2011        5887 :       gen = g;
    2012             :     }
    2013        7406 :     gel(G,3) = gen; /* replace dummy */
    2014             :   }
    2015             :   else
    2016        6160 :     gel(bid,2) = mkvec2(h,cyc);
    2017       13566 : }
    2018             : 
    2019             : static int
    2020        3472 : RgV_is_prV(GEN v)
    2021             : {
    2022        3472 :   long l = lg(v), i;
    2023        3542 :   for (i = 1; i < l; i++)
    2024        3472 :     if (!get_prid(gel(v,i))) return 0;
    2025          70 :   return 1;
    2026             : }
    2027             : 
    2028             : static int
    2029       13272 : is_nf_factor(GEN F)
    2030             : {
    2031       38605 :   return typ(F) == t_MAT && lg(F) == 3
    2032       16744 :     && RgV_is_prV(gel(F,1)) && RgV_is_ZV(gel(F,2));
    2033             : }
    2034             : 
    2035             : /* Compute [[ideal,arch], [h,[cyc],[gen]], idealfact, [liste], U]
    2036             :    flag may include nf_GEN | nf_INIT */
    2037             : GEN
    2038       13545 : Idealstar(GEN nf, GEN ideal, long flag)
    2039             : {
    2040       13545 :   pari_sp av = avma;
    2041             :   long i, j, k, nbp, R1, nbgen;
    2042       13545 :   GEN t, y, cyc, U, u1 = NULL, fa, lists, x, arch, archp, E, P, sarch, gen;
    2043             : 
    2044       13545 :   if (!nf) return znstar0(ideal, flag);
    2045       13272 :   nf = checknf(nf);
    2046       13272 :   R1 = nf_get_r1(nf);
    2047       13272 :   if (typ(ideal) == t_VEC && lg(ideal) == 3)
    2048             :   {
    2049        5285 :     arch = gel(ideal,2);
    2050        5285 :     ideal= gel(ideal,1);
    2051        5285 :     switch(typ(arch))
    2052             :     {
    2053             :       case t_VEC:
    2054        5285 :         if (lg(arch) != R1+1)
    2055           0 :           pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    2056        5285 :         archp = vec01_to_indices(arch);
    2057        5285 :         break;
    2058             :       case t_VECSMALL:
    2059           0 :         archp = arch;
    2060           0 :         arch = vec01_to_indices(archp);
    2061           0 :         break;
    2062             :       default:
    2063           0 :         pari_err_TYPE("Idealstar [incorrect archimedean component]",arch);
    2064           0 :         return NULL;
    2065             :     }
    2066        5285 :   }
    2067             :   else
    2068             :   {
    2069        7987 :     arch = zerovec(R1);
    2070        7987 :     archp = cgetg(1, t_VECSMALL);
    2071             :   }
    2072       13272 :   if (is_nf_factor(ideal))
    2073             :   {
    2074          70 :     fa = ideal;
    2075          70 :     x = idealfactorback(nf, gel(fa,1), gel(fa,2), 0);
    2076             :   }
    2077             :   else
    2078             :   {
    2079       13202 :     fa = NULL;
    2080       13202 :     x = idealhnf_shallow(nf, ideal);
    2081             :   }
    2082       13272 :   if (lg(x) == 1) pari_err_DOMAIN("Idealstar", "ideal","=",gen_0,x);
    2083       13265 :   if (typ(gcoeff(x,1,1)) != t_INT)
    2084           7 :     pari_err_DOMAIN("Idealstar","denominator(ideal)", "!=",gen_1,x);
    2085       13258 :   sarch = nfarchstar(nf, x, archp);
    2086       13258 :   if (!fa) fa = idealfactor(nf, ideal);
    2087       13258 :   P = gel(fa,1);
    2088       13258 :   E = gel(fa,2); nbp = lg(P)-1;
    2089       13258 :   if (nbp)
    2090             :   {
    2091             :     GEN h;
    2092       11683 :     long cp = 0;
    2093             :     zlog_S S;
    2094             : 
    2095       11683 :     lists = cgetg(nbp+2,t_VEC);
    2096             :     /* rough upper bound */
    2097       11683 :     nbgen = nbp + 1; for (i=1; i<=nbp; i++) nbgen += itos(gel(E,i));
    2098       11683 :     gen = cgetg(nbgen+1,t_VEC);
    2099       11683 :     nbgen = 1;
    2100       11683 :     t = (nbp==1)? NULL: x;
    2101       24598 :     for (i=1; i<=nbp; i++)
    2102             :     {
    2103       12922 :       GEN L = zprimestar(nf, gel(P,i), gel(E,i), t, archp);
    2104       12915 :       gel(lists,i) = L;
    2105       12915 :       for (j = 1; j < lg(L); j++) gel(gen, nbgen++) = gmael(L,j,3);
    2106             :     }
    2107       11676 :     gel(lists,i) = sarch;
    2108       11676 :     gel(gen, nbgen++) = gel(sarch,2); setlg(gen, nbgen);
    2109       11676 :     gen = shallowconcat1(gen); nbgen = lg(gen)-1;
    2110             : 
    2111       11676 :     h = cgetg(nbgen+1,t_MAT);
    2112       11676 :     init_zlog(&S, nbgen, P, E, archp, lists, NULL);
    2113       24591 :     for (i=1; i<=nbp; i++)
    2114             :     {
    2115       12915 :       GEN L2 = gel(lists,i);
    2116       37597 :       for (j=1; j < lg(L2); j++)
    2117             :       {
    2118       24682 :         GEN L = gel(L2,j), F = gel(L,1), G = gel(L,3);
    2119       83608 :         for (k=1; k<lg(G); k++)
    2120             :         { /* log(g^f) mod divisor */
    2121       58926 :           GEN g = gel(G,k), f = gel(F,k), a = nfpowmodideal(nf,g,f,x);
    2122      128870 :           GEN sgn = mpodd(f)? nfsign_arch(nf, g, S.archp)
    2123       69944 :                             : zero_zv(lg(S.archp)-1);
    2124       58926 :           gel(h,++cp) = ZC_neg(zlog_ind(nf, a, &S, sgn, i));
    2125       58926 :           gcoeff(h,cp,cp) = f;
    2126             :         }
    2127             :       }
    2128             :     }
    2129       16387 :     for (j=1; j<lg(archp); j++)
    2130             :     {
    2131        4711 :       gel(h,++cp) = zerocol(nbgen);
    2132        4711 :       gcoeff(h,cp,cp) = gen_2;
    2133             :     }
    2134             :     /* assert(cp == nbgen) */
    2135       11676 :     h = ZM_hnfall(h,NULL,0);
    2136       11676 :     cyc = ZM_snf_group(h, &U, (flag & nf_GEN)? &u1: NULL);
    2137             :   }
    2138             :   else
    2139             :   {
    2140        1575 :     lists = mkvec(sarch);
    2141        1575 :     gen = gel(sarch,2); nbgen = lg(gen)-1;
    2142        1575 :     cyc = const_vec(nbgen, gen_2);
    2143        1575 :     U = matid(nbgen);
    2144        1575 :     if (flag & nf_GEN) u1 = gen_1;
    2145             :   }
    2146             : 
    2147       13251 :   y = cgetg(6,t_VEC);
    2148       13251 :   gel(y,1) = mkvec2(x, arch);
    2149       13251 :   gel(y,3) = fa;
    2150       13251 :   gel(y,4) = lists;
    2151       13251 :   gel(y,5) = U;
    2152       13251 :   add_grp(nf, u1, cyc, gen, y);
    2153       13251 :   if (!(flag & nf_INIT)) y = gel(y,2);
    2154       13251 :   return gerepilecopy(av, y);
    2155             : }
    2156             : 
    2157             : /* vectors of [[cyc],[g],U.X^-1] */
    2158             : static GEN
    2159           7 : principal_units(GEN nf, GEN pr, long e)
    2160             : {
    2161           7 :   pari_sp av = avma;
    2162             :   long a;
    2163             :   GEN list, y, prb, pre;
    2164             :   ulong mask;
    2165             : 
    2166           7 :   if(DEBUGLEVEL>3) err_printf("treating pr^%ld, pr = %Ps\n",e,pr);
    2167           7 :   if (e == 1) return cgetg(1, t_VEC);
    2168           7 :   prb = idealhnf_two(nf,pr);
    2169           7 :   pre = idealpows(nf,pr,e);
    2170           7 :   list = vectrunc_init(e);
    2171           7 :   mask = quadratic_prec_mask(e);
    2172           7 :   a = 1;
    2173          42 :   while (mask > 1)
    2174             :   {
    2175          28 :     GEN pra = prb, z, U;
    2176          28 :     long b = a << 1;
    2177             : 
    2178          28 :     if (mask & 1) b--;
    2179          28 :     mask >>= 1;
    2180             :     /* compute 1 + pr^a / 1 + pr^b, 2a <= b */
    2181          28 :     if(DEBUGLEVEL>3) err_printf("  treating a = %ld, b = %ld\n",a,b);
    2182          28 :     prb = (b >= e)? pre: idealpows(nf,pr,b);
    2183          28 :     z = zidealij(pra, prb, &U);
    2184          28 :     y = mkvec3(gel(z,1), gel(z,2), U);
    2185          28 :     vectrunc_append(list, y);
    2186          28 :     a = b;
    2187             :   }
    2188           7 :   return gerepilecopy(av, list);
    2189             : }
    2190             : 
    2191             : static GEN
    2192          42 : log_prk(GEN nf, GEN a, long nbgen, GEN list, GEN prk)
    2193             : {
    2194          42 :   GEN y = zerocol(nbgen);
    2195          42 :   long i,j, iy = 1, llist = lg(list)-1;
    2196             : 
    2197         210 :   for (j = 1; j <= llist; j++)
    2198             :   {
    2199         168 :     GEN L = gel(list,j);
    2200         168 :     GEN cyc = gel(L,1), gen = gel(L,2), U = gel(L,3);
    2201         168 :     GEN e = apply_U(U, a);
    2202             :     /* here lg(e) == lg(cyc) */
    2203         420 :     for (i = 1; i < lg(cyc); i++)
    2204             :     {
    2205         252 :       GEN t = modii(negi(gel(e,i)), gel(cyc,i));
    2206         252 :       gel(y, iy++) = negi(t); if (!signe(t)) continue;
    2207          28 :       if (j != llist) a = elt_mulpow_modideal(nf, a, gel(gen,i), t, prk);
    2208             :     }
    2209             :   }
    2210          42 :   return y;
    2211             : }
    2212             : 
    2213             : /* multiplicative group (1 + pr) / (1 + pr^e) */
    2214             : GEN
    2215           7 : idealprincipalunits(GEN nf, GEN pr, long e)
    2216             : {
    2217           7 :   pari_sp av = avma;
    2218             :   long c, i, j, k, nbgen;
    2219           7 :   GEN cyc, u1 = NULL, pre, gen;
    2220             :   GEN g, EX, h, L2;
    2221           7 :   long cp = 0;
    2222             : 
    2223           7 :   nf = checknf(nf); pre = idealpows(nf, pr, e);
    2224           7 :   L2 = principal_units(nf, pr, e);
    2225           7 :   c = lg(L2); gen = cgetg(c, t_VEC);
    2226           7 :   for (j = 1; j < c; j++) gel(gen, j) = gmael(L2,j,2);
    2227           7 :   gen = shallowconcat1(gen); nbgen = lg(gen)-1;
    2228             : 
    2229           7 :   h = cgetg(nbgen+1,t_MAT);
    2230          35 :   for (j=1; j < lg(L2); j++)
    2231             :   {
    2232          28 :     GEN L = gel(L2,j), F = gel(L,1), G = gel(L,2);
    2233          70 :     for (k=1; k<lg(G); k++)
    2234             :     { /* log(g^f) mod pr^e */
    2235          42 :       GEN g = gel(G,k), f = gel(F,k), a = nfpowmodideal(nf,g,f,pre);
    2236          42 :       gel(h,++cp) = ZC_neg(log_prk(nf, a, nbgen, L2, pre));
    2237          42 :       gcoeff(h,cp,cp) = f;
    2238             :     }
    2239             :   }
    2240             :   /* assert(cp == nbgen) */
    2241           7 :   h = ZM_hnfall(h,NULL,0);
    2242           7 :   cyc = ZM_snf_group(h, NULL, &u1);
    2243           7 :   c = lg(u1); g = cgetg(c, t_VEC); EX = gel(cyc,1);
    2244          28 :   for (i=1; i<c; i++)
    2245          21 :     gel(g,i) = famat_to_nf_modideal_coprime(nf, gen, gel(u1,i), pre, EX);
    2246           7 :   return gerepilecopy(av, mkvec3(powiu(pr_norm(pr), e-1), cyc, g));
    2247             : }
    2248             : 
    2249             : /* FIXME: obsolete */
    2250             : GEN
    2251           0 : zidealstarinitgen(GEN nf, GEN ideal)
    2252           0 : { return Idealstar(nf,ideal, nf_INIT|nf_GEN); }
    2253             : GEN
    2254        4949 : zidealstarinit(GEN nf, GEN ideal)
    2255        4949 : { return Idealstar(nf,ideal, nf_INIT); }
    2256             : GEN
    2257           0 : zidealstar(GEN nf, GEN ideal)
    2258           0 : { return Idealstar(nf,ideal, nf_GEN); }
    2259             : 
    2260             : GEN
    2261         350 : idealstar0(GEN nf, GEN ideal,long flag)
    2262             : {
    2263         350 :   switch(flag)
    2264             :   {
    2265           0 :     case 0: return Idealstar(nf,ideal, nf_GEN);
    2266         315 :     case 1: return Idealstar(nf,ideal, nf_INIT);
    2267          35 :     case 2: return Idealstar(nf,ideal, nf_INIT|nf_GEN);
    2268           0 :     default: pari_err_FLAG("idealstar");
    2269             :   }
    2270           0 :   return NULL; /* not reached */
    2271             : }
    2272             : 
    2273             : void
    2274      148383 : check_nfelt(GEN x, GEN *den)
    2275             : {
    2276      148383 :   long l = lg(x), i;
    2277      148383 :   GEN t, d = NULL;
    2278      148383 :   if (typ(x) != t_COL) pari_err_TYPE("check_nfelt", x);
    2279      615305 :   for (i=1; i<l; i++)
    2280             :   {
    2281      466922 :     t = gel(x,i);
    2282      466922 :     switch (typ(t))
    2283             :     {
    2284      356379 :       case t_INT: break;
    2285             :       case t_FRAC:
    2286      110543 :         if (!d) d = gel(t,2); else d = lcmii(d, gel(t,2));
    2287      110543 :         break;
    2288           0 :       default: pari_err_TYPE("check_nfelt", x);
    2289             :     }
    2290             :   }
    2291      148383 :   *den = d;
    2292      148383 : }
    2293             : 
    2294             : GEN
    2295      413066 : vecmodii(GEN a, GEN b)
    2296             : {
    2297             :   long i, l;
    2298      413066 :   GEN c = cgetg_copy(a, &l);
    2299      413066 :   for (i = 1; i < l; i++) gel(c,i) = modii(gel(a,i), gel(b,i));
    2300      413066 :   return c;
    2301             : }
    2302             : 
    2303             : /* Given x (not necessarily integral), and bid as output by zidealstarinit,
    2304             :  * compute the vector of components on the generators bid[2].
    2305             :  * Assume (x,bid) = 1 and sgn is either NULL or nfsign_arch(x, bid) */
    2306             : GEN
    2307      199701 : ideallog_sgn(GEN nf, GEN x, GEN sgn, GEN bid)
    2308             : {
    2309             :   pari_sp av;
    2310             :   long lcyc;
    2311             :   GEN den, cyc, y;
    2312             : 
    2313      199701 :   nf = checknf(nf); checkbid(bid);
    2314      199694 :   cyc = bid_get_cyc(bid);
    2315      199694 :   lcyc = lg(cyc); if (lcyc == 1) return cgetg(1, t_COL);
    2316      199666 :   av = avma;
    2317      199666 :   if (typ(x) == t_MAT) {
    2318       42967 :     if (lg(x) == 1) return zerocol(lcyc-1); /* x = 1 */
    2319       42967 :     y = famat_zlog(nf, x, sgn, bid);
    2320       42967 :     goto END;
    2321             :   }
    2322      156699 :   x = nf_to_scalar_or_basis(nf, x);
    2323      156699 :   switch(typ(x))
    2324             :   {
    2325             :     case t_INT:
    2326        8309 :       den = NULL;
    2327        8309 :       break;
    2328             :     case t_FRAC:
    2329           7 :       den = gel(x,2);
    2330           7 :       x = gel(x,1);
    2331           7 :       break;
    2332             :     default: /* case t_COL: */
    2333      148383 :       check_nfelt(x, &den);
    2334      148383 :       if (den) x = Q_muli_to_int(x, den);
    2335             :   }
    2336      156699 :   if (den)
    2337             :   {
    2338       52377 :     x = mkmat2(mkcol2(x, den), mkcol2(gen_1, gen_m1));
    2339       52377 :     y = famat_zlog(nf, x, sgn, bid);
    2340             :   }
    2341             :   else
    2342             :   {
    2343      104322 :     zlog_S S; init_zlog_bid(&S, bid);
    2344      104322 :     y = zlog(nf, x, sgn, &S);
    2345             :   }
    2346             : END:
    2347      199659 :   y = ZM_ZC_mul(bid_get_U(bid), y);
    2348      199659 :   return gerepileupto(av, vecmodii(y, cyc));
    2349             : }
    2350             : GEN
    2351      206407 : ideallog(GEN nf, GEN x, GEN bid)
    2352             : {
    2353      206407 :   if (!nf) return Zideallog(bid, x);
    2354      199701 :   return ideallog_sgn(nf, x, NULL, bid);
    2355             : }
    2356             : 
    2357             : /*************************************************************************/
    2358             : /**                                                                     **/
    2359             : /**               JOIN BID STRUCTURES, IDEAL LISTS                      **/
    2360             : /**                                                                     **/
    2361             : /*************************************************************************/
    2362             : 
    2363             : /* bid1, bid2: for coprime modules m1 and m2 (without arch. part).
    2364             :  * Output: bid [[m1 m2,arch],[h,[cyc],[gen]],idealfact,[liste],U] for m1 m2 */
    2365             : static GEN
    2366         476 : join_bid(GEN nf, GEN bid1, GEN bid2)
    2367             : {
    2368         476 :   pari_sp av = avma;
    2369             :   long i, nbgen, lx, lx1,lx2, l1,l2;
    2370             :   GEN I1,I2, G1,G2, fa1,fa2, lists1,lists2, cyc1,cyc2;
    2371         476 :   GEN lists, fa, U, cyc, y, u1 = NULL, x, gen;
    2372             : 
    2373         476 :   I1 = bid_get_ideal(bid1);
    2374         476 :   I2 = bid_get_ideal(bid2);
    2375         476 :   if (gequal1(gcoeff(I1,1,1))) return bid2; /* frequent trivial case */
    2376         259 :   G1 = bid_get_grp(bid1);
    2377         259 :   G2 = bid_get_grp(bid2);
    2378         259 :   fa1= gel(bid1,3);
    2379         259 :   fa2= gel(bid2,3); x = idealmul(nf, I1,I2);
    2380         259 :   fa = famat_mul_shallow(fa1, fa2);
    2381         259 :   lists1 = gel(bid1,4); lx1 = lg(lists1);
    2382         259 :   lists2 = gel(bid2,4); lx2 = lg(lists2);
    2383             :   /* concat (lists1 - last elt [nfarchstar]) + lists2 */
    2384         259 :   lx = lx1+lx2-2; lists = cgetg(lx,t_VEC);
    2385         259 :   for (i=1; i<lx1-1; i++) gel(lists,i) = gel(lists1,i);
    2386         259 :   for (   ; i<lx; i++)    gel(lists,i) = gel(lists2,i-lx1+2);
    2387             : 
    2388         259 :   cyc1 = abgrp_get_cyc(G1); l1 = lg(cyc1);
    2389         259 :   cyc2 = abgrp_get_cyc(G2); l2 = lg(cyc2);
    2390         259 :   gen = (lg(G1)>3 && lg(G2)>3)? gen_1: NULL;
    2391         259 :   nbgen = l1+l2-2;
    2392         259 :   cyc = ZV_snf_group(shallowconcat(cyc1,cyc2), &U, gen? &u1: NULL);
    2393         259 :   if (nbgen) {
    2394         259 :     GEN U1 = gel(bid1,5), U2 = gel(bid2,5);
    2395         259 :     U1 = l1 == 1? zeromat(nbgen,lg(U1)-1): ZM_mul(vecslice(U, 1, l1-1),   U1);
    2396         259 :     U2 = l2 == 1? zeromat(nbgen,lg(U2)-1): ZM_mul(vecslice(U, l1, nbgen), U2);
    2397         259 :     U = shallowconcat(U1, U2);
    2398             :   }
    2399             :   else
    2400           0 :     U = zeromat(0, lx-2);
    2401             : 
    2402         259 :   if (gen)
    2403             :   {
    2404         259 :     GEN uv = zkchinese1init2(nf, I2, I1, x);
    2405         518 :     gen = shallowconcat(zkVchinese1(gel(uv,1), abgrp_get_gen(G1)),
    2406         259 :                         zkVchinese1(gel(uv,2), abgrp_get_gen(G2)));
    2407             :   }
    2408         259 :   y = cgetg(6,t_VEC);
    2409         259 :   gel(y,1) = mkvec2(x, bid_get_arch(bid1));
    2410         259 :   gel(y,3) = fa;
    2411         259 :   gel(y,4) = lists;
    2412         259 :   gel(y,5) = U;
    2413         259 :   add_grp(nf, u1, cyc, gen, y);
    2414         259 :   return gerepilecopy(av,y);
    2415             : }
    2416             : 
    2417             : typedef struct _ideal_data {
    2418             :   GEN nf, emb, L, pr, prL, arch, sgnU;
    2419             : } ideal_data;
    2420             : 
    2421             : /* z <- ( z | f(v[i])_{i=1..#v} ) */
    2422             : static void
    2423       43414 : concat_join(GEN *pz, GEN v, GEN (*f)(ideal_data*,GEN), ideal_data *data)
    2424             : {
    2425       43414 :   long i, nz, lv = lg(v);
    2426             :   GEN z, Z;
    2427       86828 :   if (lv == 1) return;
    2428       18942 :   z = *pz; nz = lg(z)-1;
    2429       18942 :   *pz = Z = cgetg(lv + nz, typ(z));
    2430       18942 :   for (i = 1; i <=nz; i++) gel(Z,i) = gel(z,i);
    2431       18942 :   Z += nz;
    2432       18942 :   for (i = 1; i < lv; i++) gel(Z,i) = f(data, gel(v,i));
    2433             : }
    2434             : static GEN
    2435         476 : join_idealinit(ideal_data *D, GEN x) {
    2436         476 :   return join_bid(D->nf, x, D->prL);
    2437             : }
    2438             : static GEN
    2439       26222 : join_ideal(ideal_data *D, GEN x) {
    2440       26222 :   return idealmulpowprime(D->nf, x, D->pr, D->L);
    2441             : }
    2442             : static GEN
    2443         455 : join_unit(ideal_data *D, GEN x) {
    2444         455 :   return mkvec2(join_idealinit(D, gel(x,1)), vconcat(gel(x,2), D->emb));
    2445             : }
    2446             : 
    2447             : /* compute matrix of zlogs of units */
    2448             : GEN
    2449          35 : zlog_units(GEN nf, GEN U, GEN sgnU, GEN bid)
    2450             : {
    2451          35 :   long j, l = lg(U);
    2452          35 :   GEN m = cgetg(l, t_MAT);
    2453          35 :   zlog_S S; init_zlog_bid(&S, bid);
    2454          98 :   for (j = 1; j < l; j++)
    2455          63 :     gel(m,j) = zlog(nf, gel(U,j), vecsmallpermute(gel(sgnU,j), S.archp), &S);
    2456          35 :   return m;
    2457             : }
    2458             : /* compute matrix of zlogs of units, assuming S.archp = [] */
    2459             : GEN
    2460         399 : zlog_units_noarch(GEN nf, GEN U, GEN bid)
    2461             : {
    2462         399 :   long j, l = lg(U);
    2463         399 :   GEN m = cgetg(l, t_MAT), empty = cgetg(1, t_COL);
    2464         399 :   zlog_S S; init_zlog_bid(&S, bid);
    2465         399 :   for (j = 1; j < l; j++) gel(m,j) = zlog(nf, gel(U,j), empty, &S);
    2466         399 :   return m;
    2467             : }
    2468             : 
    2469             : /* calcule la matrice des zlog des unites */
    2470             : static GEN
    2471          28 : zlog_unitsarch(GEN sgnU, GEN bid)
    2472             : {
    2473          28 :   GEN lists = gel(bid,4), arch = gmael(bid,1,2);
    2474          28 :   GEN listslast = gel(lists,lg(lists)-1);
    2475          28 :   GEN m = rowpermute(sgnU, vec01_to_indices(arch));
    2476          28 :   return Flm_mul(gel(listslast,3), m, 2);
    2477             : }
    2478             : 
    2479             : /*  flag & nf_GEN : generators, otherwise no
    2480             :  *  flag &2 : units, otherwise no
    2481             :  *  flag &4 : ideals in HNF, otherwise bid */
    2482             : static GEN
    2483         350 : Ideallist(GEN bnf, ulong bound, long flag)
    2484             : {
    2485         350 :   const long do_units = flag & 2, big_id = !(flag & 4);
    2486         350 :   const long istar_flag = (flag & nf_GEN) | nf_INIT;
    2487         350 :   pari_sp av, av0 = avma;
    2488             :   long i, j, l;
    2489         350 :   GEN nf, z, p, fa, id, BOUND, U, empty = cgetg(1,t_VEC);
    2490             :   forprime_t S;
    2491             :   ideal_data ID;
    2492         350 :   GEN (*join_z)(ideal_data*, GEN) =
    2493             :     do_units? &join_unit
    2494         350 :               : (big_id? &join_idealinit: &join_ideal);
    2495             : 
    2496         350 :   nf = checknf(bnf);
    2497         350 :   if ((long)bound <= 0) return empty;
    2498         350 :   id = matid(nf_get_degree(nf));
    2499         350 :   if (big_id) id = Idealstar(nf,id, istar_flag);
    2500             : 
    2501             :   /* z[i] will contain all "objects" of norm i. Depending on flag, this means
    2502             :    * an ideal, a bid, or a couple [bid, log(units)]. Such objects are stored
    2503             :    * in vectors, computed one primary component at a time; join_z
    2504             :    * reconstructs the global object */
    2505         350 :   BOUND = utoipos(bound);
    2506         350 :   z = cgetg(bound+1,t_VEC);
    2507         350 :   if (do_units) {
    2508          14 :     U = init_units(bnf);
    2509          14 :     gel(z,1) = mkvec( mkvec2(id, zlog_units_noarch(nf, U, id)) );
    2510             :   } else {
    2511         336 :     U = NULL; /* -Wall */
    2512         336 :     gel(z,1) = mkvec(id);
    2513             :   }
    2514         350 :   for (i=2; i<=(long)bound; i++) gel(z,i) = empty;
    2515         350 :   ID.nf = nf;
    2516             : 
    2517         350 :   p = cgetipos(3);
    2518         350 :   u_forprime_init(&S, 2, bound);
    2519         350 :   av = avma;
    2520        5726 :   while ((p[2] = u_forprime_next(&S)))
    2521             :   {
    2522        5026 :     if (DEBUGLEVEL>1) { err_printf("%ld ",p[2]); err_flush(); }
    2523        5026 :     fa = idealprimedec_limit_norm(nf, p, BOUND);
    2524       10073 :     for (j=1; j<lg(fa); j++)
    2525             :     {
    2526        5047 :       GEN pr = gel(fa,j), z2;
    2527        5047 :       ulong q, Q = upowuu(p[2], pr_get_f(pr));
    2528             : 
    2529        5047 :       z2 = leafcopy(z);
    2530        5047 :       q = Q;
    2531        5047 :       ID.pr = ID.prL = pr;
    2532       13370 :       for (l=1; Q <= bound; l++, Q *= q) /* add pr^l */
    2533             :       {
    2534             :         ulong iQ;
    2535        8323 :         ID.L = utoipos(l);
    2536        8323 :         if (big_id) {
    2537         217 :           if (l > 1) ID.prL = idealpow(nf,pr,ID.L);
    2538         217 :           ID.prL = Idealstar(nf,ID.prL, istar_flag);
    2539         217 :           if (do_units) ID.emb = zlog_units_noarch(nf, U, ID.prL);
    2540             :         }
    2541       51737 :         for (iQ = Q,i = 1; iQ <= bound; iQ += Q,i++)
    2542       43414 :           concat_join(&gel(z,iQ), gel(z2,i), join_z, &ID);
    2543             :       }
    2544             :     }
    2545        5026 :     if (gc_needed(av,1))
    2546             :     {
    2547           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"Ideallist");
    2548           0 :       z = gerepilecopy(av, z);
    2549             :     }
    2550             :   }
    2551         763 :   if (do_units) for (i = 1; i < lg(z); i++)
    2552             :   {
    2553         413 :     GEN s = gel(z,i);
    2554         413 :     long l = lg(s);
    2555         882 :     for (j = 1; j < l; j++) {
    2556         469 :       GEN v = gel(s,j), bid = gel(v,1);
    2557         469 :       gel(v,2) = ZM_mul(bid_get_U(bid), gel(v,2));
    2558             :     }
    2559             :   }
    2560         350 :   return gerepilecopy(av0, z);
    2561             : }
    2562             : GEN
    2563          28 : ideallist0(GEN bnf,long bound, long flag) {
    2564          28 :   if (flag<0 || flag>4) pari_err_FLAG("ideallist");
    2565          28 :   return Ideallist(bnf,bound,flag);
    2566             : }
    2567             : GEN
    2568         322 : ideallist(GEN bnf,long bound) { return Ideallist(bnf,bound,4); }
    2569             : 
    2570             : /* bid1 = for module m1 (without arch. part), arch = archimedean part.
    2571             :  * Output: bid [[m1,arch],[h,[cyc],[gen]],idealfact,[liste],U] for m1.arch */
    2572             : static GEN
    2573          56 : join_bid_arch(GEN nf, GEN bid1, GEN arch)
    2574             : {
    2575          56 :   pari_sp av = avma;
    2576             :   GEN f1, G1, fa1, U;
    2577          56 :   GEN lists, cyc, y, u1 = NULL, x, sarch, gen;
    2578             : 
    2579          56 :   checkbid(bid1);
    2580          56 :   f1 = gel(bid1,1); G1 = gel(bid1,2); fa1 = gel(bid1,3);
    2581          56 :   x = gel(f1,1);
    2582          56 :   sarch = nfarchstar(nf, x, arch);
    2583          56 :   lists = leafcopy(gel(bid1,4));
    2584          56 :   gel(lists,lg(lists)-1) = sarch;
    2585             : 
    2586          56 :   gen = (lg(G1)>3)? gen_1: NULL;
    2587          56 :   cyc = diagonal_shallow(shallowconcat(gel(G1,2), gel(sarch,1)));
    2588          56 :   cyc = ZM_snf_group(cyc, &U, gen? &u1: NULL);
    2589          56 :   if (gen) gen = shallowconcat(gel(G1,3), gel(sarch,2));
    2590          56 :   y = cgetg(6,t_VEC);
    2591          56 :   gel(y,1) = mkvec2(x, arch);
    2592          56 :   gel(y,3) = fa1;
    2593          56 :   gel(y,4) = lists;
    2594          56 :   gel(y,5) = U;
    2595          56 :   add_grp(nf, u1, cyc, gen, y);
    2596          56 :   return gerepilecopy(av,y);
    2597             : }
    2598             : static GEN
    2599          56 : join_arch(ideal_data *D, GEN x) {
    2600          56 :   return join_bid_arch(D->nf, x, D->arch);
    2601             : }
    2602             : static GEN
    2603          28 : join_archunit(ideal_data *D, GEN x) {
    2604          28 :   GEN bid = join_arch(D, gel(x,1)), u = gel(x,2);
    2605          28 :   u = ZM_mul(bid_get_U(bid), vconcat(u, zm_to_ZM(zlog_unitsarch(D->sgnU,bid))));
    2606          28 :   return mkvec2(bid, u);
    2607             : }
    2608             : 
    2609             : /* L from ideallist, add archimedean part */
    2610             : GEN
    2611          14 : ideallistarch(GEN bnf, GEN L, GEN arch)
    2612             : {
    2613             :   pari_sp av;
    2614          14 :   long i, j, l = lg(L), lz;
    2615             :   GEN v, z, V;
    2616             :   ideal_data ID;
    2617             :   GEN (*join_z)(ideal_data*, GEN);
    2618             : 
    2619          14 :   if (typ(L) != t_VEC) pari_err_TYPE("ideallistarch",L);
    2620          14 :   if (l == 1) return cgetg(1,t_VEC);
    2621          14 :   z = gel(L,1);
    2622          14 :   if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    2623          14 :   z = gel(z,1); /* either a bid or [bid,U] */
    2624          14 :   if (lg(z) == 3) { /* the latter: do units */
    2625           7 :     if (typ(z) != t_VEC) pari_err_TYPE("ideallistarch",z);
    2626           7 :     ID.sgnU = nfsign_units(bnf, NULL, 1);
    2627           7 :     join_z = &join_archunit;
    2628             :   } else
    2629           7 :     join_z = &join_arch;
    2630          14 :   ID.nf = checknf(bnf);
    2631          14 :   ID.arch = vec01_to_indices(arch);
    2632          14 :   av = avma; V = cgetg(l, t_VEC);
    2633          70 :   for (i = 1; i < l; i++)
    2634             :   {
    2635          56 :     z = gel(L,i); lz = lg(z);
    2636          56 :     gel(V,i) = v = cgetg(lz,t_VEC);
    2637          56 :     for (j=1; j<lz; j++) gel(v,j) = join_z(&ID, gel(z,j));
    2638             :   }
    2639          14 :   return gerepilecopy(av,V);
    2640             : }

Generated by: LCOV version 1.11