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 19350-bd5f220) Lines: 1382 1472 93.9 %
Date: 2016-08-24 06:11:24 Functions: 136 143 95.1 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.11