Code coverage tests

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

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

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

LCOV - code coverage report
Current view: top level - basemath - RgX.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 25406-bf255ab81b) Lines: 1531 1706 89.7 %
Date: 2020-06-04 05:59:24 Functions: 181 197 91.9 %
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             : #include "pari.h"
      15             : #include "paripriv.h"
      16             : 
      17             : /*******************************************************************/
      18             : /*                                                                 */
      19             : /*                         GENERIC                                 */
      20             : /*                                                                 */
      21             : /*******************************************************************/
      22             : 
      23             : /* Return optimal parameter l for the evaluation of n/m polynomials of degree d
      24             :    Fractional values can be used if the evaluations are done with different
      25             :    accuracies, and thus have different weights.
      26             :  */
      27             : long
      28     4154607 : brent_kung_optpow(long d, long n, long m)
      29             : {
      30             :   long p, r;
      31     4154607 :   long pold=1, rold=n*(d-1);
      32    25454575 :   for(p=2; p<=d; p++)
      33             :   {
      34    21299968 :     r = m*(p-1) + n*((d-1)/p);
      35    21299968 :     if (r<rold) { pold=p; rold=r; }
      36             :   }
      37     4154607 :   return pold;
      38             : }
      39             : 
      40             : static GEN
      41     7098668 : gen_RgXQ_eval_powers(GEN P, GEN V, long a, long n, void *E, const struct bb_algebra *ff,
      42             :                                            GEN cmul(void *E, GEN P, long a, GEN x))
      43             : {
      44     7098668 :   pari_sp av = avma;
      45             :   long i;
      46     7098668 :   GEN z = cmul(E,P,a,ff->one(E));
      47     7098231 :   if (!z) z = gen_0;
      48    50917367 :   for (i=1; i<=n; i++)
      49             :   {
      50    43819460 :     GEN t = cmul(E,P,a+i,gel(V,i+1));
      51    43822866 :     if (t) {
      52    32601688 :       z = ff->add(E, z, t);
      53    32581553 :       if (gc_needed(av,2)) z = gerepileupto(av, z);
      54             :     }
      55             :   }
      56     7097907 :   return ff->red(E,z);
      57             : }
      58             : 
      59             : /* Brent & Kung
      60             :  * (Fast algorithms for manipulating formal power series, JACM 25:581-595, 1978)
      61             :  *
      62             :  * V as output by FpXQ_powers(x,l,T,p). For optimal performance, l is as given
      63             :  * by brent_kung_optpow */
      64             : GEN
      65     5899932 : gen_bkeval_powers(GEN P, long d, GEN V, void *E, const struct bb_algebra *ff,
      66             :                                      GEN cmul(void *E, GEN P, long a, GEN x))
      67             : {
      68     5899932 :   pari_sp av = avma;
      69     5899932 :   long l = lg(V)-1;
      70             :   GEN z, u;
      71             : 
      72     5899932 :   if (d < 0) return ff->zero(E);
      73     5581985 :   if (d < l) return gerepileupto(av, gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul));
      74      651865 :   if (l<2) pari_err_DOMAIN("gen_RgX_bkeval_powers", "#powers", "<",gen_2,V);
      75      651865 :   if (DEBUGLEVEL>=8)
      76             :   {
      77           0 :     long cnt = 1 + (d - l) / (l-1);
      78           0 :     err_printf("RgX_RgXQV_eval(%ld/%ld): %ld RgXQ_mul\n", d, l-1, cnt);
      79             :   }
      80      651865 :   d -= l;
      81      651865 :   z = gen_RgXQ_eval_powers(P,V,d+1,l-1,E,ff,cmul);
      82     1517058 :   while (d >= l-1)
      83             :   {
      84      865215 :     d -= l-1;
      85      865215 :     u = gen_RgXQ_eval_powers(P,V,d+1,l-2,E,ff,cmul);
      86      865324 :     z = ff->add(E,u, ff->mul(E,z,gel(V,l)));
      87      865177 :     if (gc_needed(av,2))
      88          91 :       z = gerepileupto(av, z);
      89             :   }
      90      651843 :   u = gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul);
      91      651850 :   z = ff->add(E,u, ff->mul(E,z,gel(V,d+2)));
      92      651849 :   return gerepileupto(av, ff->red(E,z));
      93             : }
      94             : 
      95             : GEN
      96      424875 : gen_bkeval(GEN Q, long d, GEN x, int use_sqr, void *E, const struct bb_algebra *ff,
      97             :                                       GEN cmul(void *E, GEN P, long a, GEN x))
      98             : {
      99      424875 :   pari_sp av = avma;
     100             :   GEN z, V;
     101             :   long rtd;
     102      424875 :   if (d < 0) return ff->zero(E);
     103      424308 :   rtd = (long) sqrt((double)d);
     104      424308 :   V = gen_powers(x,rtd,use_sqr,E,ff->sqr,ff->mul,ff->one);
     105      424314 :   z = gen_bkeval_powers(Q, d, V, E, ff, cmul);
     106      424309 :   return gerepileupto(av, z);
     107             : }
     108             : 
     109             : static GEN
     110     1092503 : _gen_nored(void *E, GEN x) { (void)E; return x; }
     111             : static GEN
     112    13974287 : _gen_add(void *E, GEN x, GEN y) { (void)E; return gadd(x, y); }
     113             : static GEN
     114           0 : _gen_sub(void *E, GEN x, GEN y) { (void)E; return gsub(x, y); }
     115             : static GEN
     116      941430 : _gen_mul(void *E, GEN x, GEN y) { (void)E; return gmul(x, y); }
     117             : static GEN
     118      308270 : _gen_sqr(void *E, GEN x) { (void)E; return gsqr(x); }
     119             : static GEN
     120     1114947 : _gen_one(void *E) { (void)E; return gen_1; }
     121             : static GEN
     122       21140 : _gen_zero(void *E) { (void)E; return gen_0; }
     123             : 
     124             : static struct bb_algebra Rg_algebra = { _gen_nored, _gen_add, _gen_sub,
     125             :               _gen_mul, _gen_sqr,_gen_one,_gen_zero };
     126             : 
     127             : static GEN
     128      479461 : _gen_cmul(void *E, GEN P, long a, GEN x)
     129      479461 : {(void)E; return gmul(gel(P,a+2), x);}
     130             : 
     131             : GEN
     132      155400 : RgX_RgV_eval(GEN Q, GEN x)
     133             : {
     134      155400 :   return gen_bkeval_powers(Q, degpol(Q), x, NULL, &Rg_algebra, _gen_cmul);
     135             : }
     136             : 
     137             : GEN
     138           0 : RgX_Rg_eval_bk(GEN Q, GEN x)
     139             : {
     140           0 :   return gen_bkeval(Q, degpol(Q), x, 1, NULL, &Rg_algebra, _gen_cmul);
     141             : }
     142             : 
     143             : GEN
     144        2912 : RgXV_RgV_eval(GEN Q, GEN x)
     145             : {
     146        2912 :   long i, l = lg(Q), vQ = gvar(Q);
     147        2912 :   GEN v = cgetg(l, t_VEC);
     148      243327 :   for (i = 1; i < l; i++)
     149             :   {
     150      240415 :     GEN Qi = gel(Q, i);
     151      240415 :     gel(v, i) = typ(Qi)==t_POL && varn(Qi)==vQ? RgX_RgV_eval(Qi, x): gcopy(Qi);
     152             :   }
     153        2912 :   return v;
     154             : }
     155             : 
     156             : GEN
     157        5040 : RgX_homogenous_evalpow(GEN P, GEN A, GEN B)
     158             : {
     159        5040 :   pari_sp av = avma;
     160             :   long d, i, v;
     161             :   GEN s;
     162        5040 :   if (typ(P)!=t_POL)
     163        1547 :     return mkvec2(P, gen_1);
     164        3493 :   d = degpol(P); v = varn(A);
     165        3493 :   s = scalarpol_shallow(gel(P, d+2), v);
     166       17955 :   for (i = d-1; i >= 0; i--)
     167             :   {
     168       14462 :     s = gadd(gmul(s, A), gmul(gel(B,d+1-i), gel(P,i+2)));
     169       14462 :     if (gc_needed(av,1))
     170             :     {
     171           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_homogenous_eval(%ld)",i);
     172           0 :       s = gerepileupto(av, s);
     173             :     }
     174             :   }
     175        3493 :   s = gerepileupto(av, s);
     176        3493 :   return mkvec2(s, gel(B,d+1));
     177             : }
     178             : 
     179             : GEN
     180        1792 : RgXQX_homogenous_evalpow(GEN P, GEN A, GEN B, GEN T)
     181             : {
     182        1792 :   pari_sp av = avma;
     183        1792 :   long i, d = degpol(P), v = varn(A);
     184             :   GEN s;
     185        1792 :   if (signe(P)==0) return mkvec2(pol_0(v), pol_1(v));
     186        1568 :   s = scalarpol_shallow(gel(P, d+2), v);
     187        5229 :   for (i = d-1; i >= 0; i--)
     188             :   {
     189        3661 :     s = RgX_add(RgXQX_mul(s, A, T), RgXQX_RgXQ_mul(gel(B,d+1-i), gel(P,i+2), T));
     190        3661 :     if (gc_needed(av,1))
     191             :     {
     192           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_homogenous_eval(%ld)",i);
     193           0 :       s = gerepileupto(av, s);
     194             :     }
     195             :   }
     196        1568 :   s = gerepileupto(av, s);
     197        1568 :   return mkvec2(s, gel(B,d+1));
     198             : }
     199             : 
     200             : const struct bb_algebra *
     201      142984 : get_Rg_algebra(void)
     202             : {
     203      142984 :   return &Rg_algebra;
     204             : }
     205             : 
     206             : static struct bb_ring Rg_ring = {  _gen_add, _gen_mul, _gen_sqr };
     207             : 
     208             : static GEN
     209       17780 : _RgX_divrem(void *E, GEN x, GEN y, GEN *r)
     210             : {
     211             :   (void) E;
     212       17780 :   return RgX_divrem(x, y, r);
     213             : }
     214             : 
     215             : GEN
     216        4382 : RgX_digits(GEN x, GEN T)
     217             : {
     218        4382 :   pari_sp av = avma;
     219        4382 :   long d = degpol(T), n = (lgpol(x)+d-1)/d;
     220        4382 :   GEN z = gen_digits(x,T,n,NULL, &Rg_ring, _RgX_divrem);
     221        4382 :   return gerepileupto(av, z);
     222             : }
     223             : 
     224             : /*******************************************************************/
     225             : /*                                                                 */
     226             : /*                         RgX                                     */
     227             : /*                                                                 */
     228             : /*******************************************************************/
     229             : 
     230             : long
     231    22202580 : RgX_equal(GEN x, GEN y)
     232             : {
     233    22202580 :   long i = lg(x);
     234             : 
     235    22202580 :   if (i != lg(y)) return 0;
     236    96396552 :   for (i--; i > 1; i--)
     237    74323346 :     if (!gequal(gel(x,i),gel(y,i))) return 0;
     238    22073206 :   return 1;
     239             : }
     240             : 
     241             : /* Returns 1 in the base ring over which x is defined */
     242             : /* HACK: this also works for t_SER */
     243             : GEN
     244     3468977 : Rg_get_1(GEN x)
     245             : {
     246             :   GEN p, T;
     247     3468977 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     248     3468977 :   if (RgX_type_is_composite(tx))
     249       34048 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     250     3468977 :   switch(tx)
     251             :   {
     252         357 :     case t_INTMOD: retmkintmod(is_pm1(p)? gen_0: gen_1, icopy(p));
     253           7 :     case t_PADIC: return cvtop(gen_1, p, lx);
     254         735 :     case t_FFELT: return FF_1(T);
     255     3467878 :     default: return gen_1;
     256             :   }
     257             : }
     258             : /* Returns 0 in the base ring over which x is defined */
     259             : /* HACK: this also works for t_SER */
     260             : GEN
     261     6972315 : Rg_get_0(GEN x)
     262             : {
     263             :   GEN p, T;
     264     6972315 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     265     6972315 :   if (RgX_type_is_composite(tx))
     266       51310 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     267     6972315 :   switch(tx)
     268             :   {
     269         441 :     case t_INTMOD: retmkintmod(gen_0, icopy(p));
     270           0 :     case t_PADIC: return zeropadic(p, lx);
     271         266 :     case t_FFELT: return FF_zero(T);
     272     6971608 :     default: return gen_0;
     273             :   }
     274             : }
     275             : 
     276             : GEN
     277        4298 : QX_ZXQV_eval(GEN P, GEN V, GEN dV)
     278             : {
     279        4298 :   long i, n = degpol(P);
     280             :   GEN z, dz, dP;
     281        4298 :   if (n < 0) return gen_0;
     282        4298 :   P = Q_remove_denom(P, &dP);
     283        4298 :   z = gel(P,2); if (n == 0) return icopy(z);
     284        2345 :   if (dV) z = mulii(dV, z); /* V[1] = dV */
     285        2345 :   z = ZX_Z_add_shallow(ZX_Z_mul(gel(V,2),gel(P,3)), z);
     286        4949 :   for (i=2; i<=n; i++) z = ZX_add(ZX_Z_mul(gel(V,i+1),gel(P,2+i)), z);
     287        2345 :   dz = mul_denom(dP, dV);
     288        2345 :   return dz? RgX_Rg_div(z, dz): z;
     289             : }
     290             : 
     291             : /* Return P(h * x), not memory clean */
     292             : GEN
     293        6930 : RgX_unscale(GEN P, GEN h)
     294             : {
     295        6930 :   long i, l = lg(P);
     296        6930 :   GEN hi = gen_1, Q = cgetg(l, t_POL);
     297        6930 :   Q[1] = P[1];
     298        6930 :   if (l == 2) return Q;
     299        6902 :   gel(Q,2) = gcopy(gel(P,2));
     300       37184 :   for (i=3; i<l; i++)
     301             :   {
     302       30282 :     hi = gmul(hi,h);
     303       30282 :     gel(Q,i) = gmul(gel(P,i), hi);
     304             :   }
     305        6902 :   return Q;
     306             : }
     307             : /* P a ZX, Return P(h * x), not memory clean; optimize for h = -1 */
     308             : GEN
     309      372913 : ZX_z_unscale(GEN P, long h)
     310             : {
     311      372913 :   long i, l = lg(P);
     312      372913 :   GEN Q = cgetg(l, t_POL);
     313      372913 :   Q[1] = P[1];
     314      372913 :   if (l == 2) return Q;
     315      372913 :   gel(Q,2) = gel(P,2);
     316      372913 :   if (l == 3) return Q;
     317      372913 :   if (h == -1)
     318      323288 :     for (i = 3; i < l; i++)
     319             :     {
     320      311623 :       gel(Q,i) = negi(gel(P,i));
     321      311623 :       if (++i == l) break;
     322      307844 :       gel(Q,i) = gel(P,i);
     323             :     }
     324             :   else
     325             :   {
     326             :     GEN hi;
     327      357469 :     gel(Q,3) = mulis(gel(P,3), h);
     328      357469 :     hi = sqrs(h);
     329     1724016 :     for (i = 4; i < l; i++)
     330             :     {
     331     1366547 :       gel(Q,i) = mulii(gel(P,i), hi);
     332     1366547 :       if (i != l-1) hi = mulis(hi,h);
     333             :     }
     334             :   }
     335      372913 :   return Q;
     336             : }
     337             : /* P a ZX, h a t_INT. Return P(h * x), not memory clean; optimize for h = -1 */
     338             : GEN
     339        7273 : ZX_unscale(GEN P, GEN h)
     340             : {
     341             :   long i, l;
     342             :   GEN Q, hi;
     343        7273 :   i = itos_or_0(h); if (i) return ZX_z_unscale(P, i);
     344           7 :   l = lg(P); Q = cgetg(l, t_POL);
     345           7 :   Q[1] = P[1];
     346           7 :   if (l == 2) return Q;
     347           7 :   gel(Q,2) = gel(P,2);
     348           7 :   if (l == 3) return Q;
     349           7 :   hi = h;
     350           7 :   gel(Q,3) = mulii(gel(P,3), hi);
     351          70 :   for (i = 4; i < l; i++)
     352             :   {
     353          63 :     hi = mulii(hi,h);
     354          63 :     gel(Q,i) = mulii(gel(P,i), hi);
     355             :   }
     356           7 :   return Q;
     357             : }
     358             : /* P a ZX. Return P(x << n), not memory clean */
     359             : GEN
     360       64007 : ZX_unscale2n(GEN P, long n)
     361             : {
     362       64007 :   long i, ni = n, l = lg(P);
     363       64007 :   GEN Q = cgetg(l, t_POL);
     364       64007 :   Q[1] = P[1];
     365       64007 :   if (l == 2) return Q;
     366       64007 :   gel(Q,2) = gel(P,2);
     367       64007 :   if (l == 3) return Q;
     368       64007 :   gel(Q,3) = shifti(gel(P,3), ni);
     369      254666 :   for (i=4; i<l; i++)
     370             :   {
     371      190659 :     ni += n;
     372      190659 :     gel(Q,i) = shifti(gel(P,i), ni);
     373             :   }
     374       64007 :   return Q;
     375             : }
     376             : /* P(h*X) / h, assuming h | P(0), i.e. the result is a ZX */
     377             : GEN
     378        1190 : ZX_unscale_div(GEN P, GEN h)
     379             : {
     380        1190 :   long i, l = lg(P);
     381        1190 :   GEN hi, Q = cgetg(l, t_POL);
     382        1190 :   Q[1] = P[1];
     383        1190 :   if (l == 2) return Q;
     384        1190 :   gel(Q,2) = diviiexact(gel(P,2), h);
     385        1190 :   if (l == 3) return Q;
     386        1190 :   gel(Q,3) = gel(P,3);
     387        1190 :   if (l == 4) return Q;
     388        1190 :   hi = h;
     389        1190 :   gel(Q,4) = mulii(gel(P,4), hi);
     390        5257 :   for (i=5; i<l; i++)
     391             :   {
     392        4067 :     hi = mulii(hi,h);
     393        4067 :     gel(Q,i) = mulii(gel(P,i), hi);
     394             :   }
     395        1190 :   return Q;
     396             : }
     397             : 
     398             : GEN
     399         658 : RgXV_unscale(GEN v, GEN h)
     400             : {
     401             :   long i, l;
     402             :   GEN w;
     403         658 :   if (!h || isint1(h)) return v;
     404         350 :   w = cgetg_copy(v, &l);
     405        1435 :   for (i=1; i<l; i++) gel(w,i) = RgX_unscale(gel(v,i), h);
     406         350 :   return w;
     407             : }
     408             : 
     409             : /* Return h^degpol(P) P(x / h), not memory clean */
     410             : GEN
     411      488152 : RgX_rescale(GEN P, GEN h)
     412             : {
     413      488152 :   long i, l = lg(P);
     414      488152 :   GEN Q = cgetg(l,t_POL), hi = h;
     415      488152 :   gel(Q,l-1) = gel(P,l-1);
     416     1767799 :   for (i=l-2; i>=2; i--)
     417             :   {
     418     1765510 :     gel(Q,i) = gmul(gel(P,i), hi);
     419     1765510 :     if (i == 2) break;
     420     1279647 :     hi = gmul(hi,h);
     421             :   }
     422      488152 :   Q[1] = P[1]; return Q;
     423             : }
     424             : 
     425             : /* A(X^d) --> A(X) */
     426             : GEN
     427      158016 : RgX_deflate(GEN x0, long d)
     428             : {
     429             :   GEN z, y, x;
     430      158016 :   long i,id, dy, dx = degpol(x0);
     431      158016 :   if (d == 1 || dx <= 0) return leafcopy(x0);
     432       83696 :   dy = dx/d;
     433       83696 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     434       83696 :   z = y + 2;
     435       83696 :   x = x0+ 2;
     436      371156 :   for (i=id=0; i<=dy; i++,id+=d) gel(z,i) = gel(x,id);
     437       83696 :   return y;
     438             : }
     439             : 
     440             : /* F a t_RFRAC */
     441             : long
     442          70 : rfrac_deflate_order(GEN F)
     443             : {
     444          70 :   GEN N = gel(F,1), D = gel(F,2);
     445          70 :   long m = (degpol(D) <= 0)? 0: RgX_deflate_order(D);
     446          70 :   if (m == 1) return 1;
     447          28 :   if (typ(N) == t_POL && varn(N) == varn(D))
     448          14 :     m = cgcd(m, RgX_deflate_order(N));
     449          28 :   return m;
     450             : }
     451             : /* F a t_RFRAC */
     452             : GEN
     453          70 : rfrac_deflate_max(GEN F, long *m)
     454             : {
     455          70 :   *m = rfrac_deflate_order(F);
     456          70 :   return rfrac_deflate(F, *m);
     457             : }
     458             : /* F a t_RFRAC */
     459             : GEN
     460          70 : rfrac_deflate(GEN F, long m)
     461             : {
     462          70 :   GEN N = gel(F,1), D = gel(F,2);
     463          70 :   if (m == 1) return F;
     464          28 :   if (typ(N) == t_POL && varn(N) == varn(D)) N = RgX_deflate(N, m);
     465          28 :   D = RgX_deflate(D, m); return mkrfrac(N, D);
     466             : }
     467             : 
     468             : /* return x0(X^d) */
     469             : GEN
     470      601268 : RgX_inflate(GEN x0, long d)
     471             : {
     472      601268 :   long i, id, dy, dx = degpol(x0);
     473      601268 :   GEN x = x0 + 2, z, y;
     474      601268 :   if (dx <= 0) return leafcopy(x0);
     475      593603 :   dy = dx*d;
     476      593603 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     477      593603 :   z = y + 2;
     478    21657842 :   for (i=0; i<=dy; i++) gel(z,i) = gen_0;
     479     9727312 :   for (i=id=0; i<=dx; i++,id+=d) gel(z,id) = gel(x,i);
     480      593603 :   return y;
     481             : }
     482             : 
     483             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     484             : static GEN
     485     1698742 : RgX_translate_basecase(GEN P, GEN c)
     486             : {
     487     1698742 :   pari_sp av = avma;
     488             :   GEN Q, R;
     489             :   long i, k, n;
     490             : 
     491     1698742 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     492     1696338 :   Q = leafcopy(P);
     493     1696338 :   R = Q+2; n = degpol(P);
     494     1696338 :   if (isint1(c))
     495             :   {
     496        9023 :     for (i=1; i<=n; i++)
     497             :     {
     498       29505 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gel(R,k+1));
     499        6657 :       if (gc_needed(av,2))
     500             :       {
     501           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(1), i = %ld/%ld", i,n);
     502           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     503             :       }
     504             :     }
     505             :   }
     506     1693972 :   else if (isintm1(c))
     507             :   {
     508      162860 :     for (i=1; i<=n; i++)
     509             :     {
     510     1004649 :       for (k=n-i; k<n; k++) gel(R,k) = gsub(gel(R,k), gel(R,k+1));
     511      138531 :       if (gc_needed(av,2))
     512             :       {
     513           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(-1), i = %ld/%ld", i,n);
     514           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     515             :       }
     516             :     }
     517             :   }
     518             :   else
     519             :   {
     520     5552145 :     for (i=1; i<=n; i++)
     521             :     {
     522    12832008 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gmul(c, gel(R,k+1)));
     523     3882502 :       if (gc_needed(av,2))
     524             :       {
     525           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate, i = %ld/%ld", i,n);
     526           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     527             :       }
     528             :     }
     529             :   }
     530     1696338 :   return gerepilecopy(av, Q);
     531             : }
     532             : GEN
     533     1699036 : RgX_translate(GEN P, GEN c)
     534             : {
     535     1699036 :   pari_sp av = avma;
     536     1699036 :   long n = degpol(P);
     537     1699036 :   if (n < 40)
     538     1698742 :     return RgX_translate_basecase(P, c);
     539             :   else
     540             :   {
     541         294 :     long d = n >> 1;
     542         294 :     GEN Q = RgX_translate(RgX_shift_shallow(P, -d), c);
     543         294 :     GEN R = RgX_translate(RgXn_red_shallow(P, d), c);
     544         294 :     GEN S = gpowgs(deg1pol_shallow(gen_1, c, varn(P)), d);
     545         294 :     return gerepileupto(av, RgX_add(RgX_mul(Q, S), R));
     546             :   }
     547             : }
     548             : 
     549             : /* return lift( P(X + c) ) using Horner, c in R[y]/(T) */
     550             : GEN
     551       21280 : RgXQX_translate(GEN P, GEN c, GEN T)
     552             : {
     553       21280 :   pari_sp av = avma;
     554             :   GEN Q, R;
     555             :   long i, k, n;
     556             : 
     557       21280 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     558       20951 :   Q = leafcopy(P);
     559       20951 :   R = Q+2; n = degpol(P);
     560       72842 :   for (i=1; i<=n; i++)
     561             :   {
     562      218806 :     for (k=n-i; k<n; k++)
     563             :     {
     564      166915 :       pari_sp av2 = avma;
     565      166915 :       gel(R,k) = gerepileupto(av2,
     566      166915 :                    RgX_rem(gadd(gel(R,k), gmul(c, gel(R,k+1))), T));
     567             :     }
     568       51891 :     if (gc_needed(av,2))
     569             :     {
     570           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXQX_translate, i = %ld/%ld", i,n);
     571           0 :       Q = gerepilecopy(av, Q); R = Q+2;
     572             :     }
     573             :   }
     574       20951 :   return gerepilecopy(av, Q);
     575             : }
     576             : 
     577             : /********************************************************************/
     578             : /**                                                                **/
     579             : /**                          CONVERSIONS                           **/
     580             : /**                       (not memory clean)                       **/
     581             : /**                                                                **/
     582             : /********************************************************************/
     583             : /* to INT / FRAC / (POLMOD mod T), not memory clean because T not copied,
     584             :  * but everything else is */
     585             : static GEN
     586       39466 : QXQ_to_mod(GEN x, GEN T)
     587             : {
     588             :   long d;
     589       39466 :   switch(typ(x))
     590             :   {
     591       16674 :     case t_INT:  return icopy(x);
     592         805 :     case t_FRAC: return gcopy(x);
     593       21987 :     case t_POL:
     594       21987 :       d = degpol(x);
     595       21987 :       if (d < 0) return gen_0;
     596       21182 :       if (d == 0) return gcopy(gel(x,2));
     597       20552 :       return mkpolmod(RgX_copy(x), T);
     598           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     599             :              return NULL;/* LCOV_EXCL_LINE */
     600             :   }
     601             : }
     602             : /* pure shallow version */
     603             : GEN
     604     1216384 : QXQ_to_mod_shallow(GEN x, GEN T)
     605             : {
     606             :   long d;
     607     1216384 :   switch(typ(x))
     608             :   {
     609      459529 :     case t_INT:
     610      459529 :     case t_FRAC: return x;
     611      756855 :     case t_POL:
     612      756855 :       d = degpol(x);
     613      756855 :       if (d < 0) return gen_0;
     614      579480 :       if (d == 0) return gel(x,2);
     615      556869 :       return mkpolmod(x, T);
     616           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     617             :              return NULL;/* LCOV_EXCL_LINE */
     618             :   }
     619             : }
     620             : /* T a ZX, z lifted from (Q[Y]/(T(Y)))[X], apply QXQ_to_mod to all coeffs.
     621             :  * Not memory clean because T not copied, but everything else is */
     622             : static GEN
     623        7994 : QXQX_to_mod(GEN z, GEN T)
     624             : {
     625        7994 :   long i,l = lg(z);
     626        7994 :   GEN x = cgetg(l,t_POL);
     627       40873 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod(gel(z,i), T);
     628        7994 :   x[1] = z[1]; return normalizepol_lg(x,l);
     629             : }
     630             : /* pure shallow version */
     631             : GEN
     632      188930 : QXQX_to_mod_shallow(GEN z, GEN T)
     633             : {
     634      188930 :   long i,l = lg(z);
     635      188930 :   GEN x = cgetg(l,t_POL);
     636      855009 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod_shallow(gel(z,i), T);
     637      188930 :   x[1] = z[1]; return normalizepol_lg(x,l);
     638             : }
     639             : /* Apply QXQX_to_mod to all entries. Memory-clean ! */
     640             : GEN
     641        2527 : QXQXV_to_mod(GEN V, GEN T)
     642             : {
     643        2527 :   long i, l = lg(V);
     644        2527 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     645       10521 :   for (i=1;i<l; i++) gel(z,i) = QXQX_to_mod(gel(V,i), T);
     646        2527 :   return z;
     647             : }
     648             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     649             : GEN
     650        5124 : QXQV_to_mod(GEN V, GEN T)
     651             : {
     652        5124 :   long i, l = lg(V);
     653        5124 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     654       11711 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod(gel(V,i), T);
     655        5124 :   return z;
     656             : }
     657             : 
     658             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     659             : GEN
     660       24465 : QXQC_to_mod_shallow(GEN V, GEN T)
     661             : {
     662       24465 :   long i, l = lg(V);
     663       24465 :   GEN z = cgetg(l, t_COL);
     664      574770 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod_shallow(gel(V,i), T);
     665       24465 :   return z;
     666             : }
     667             : 
     668             : GEN
     669        7476 : QXQM_to_mod_shallow(GEN V, GEN T)
     670             : {
     671        7476 :   long i, l = lg(V);
     672        7476 :   GEN z = cgetg(l, t_MAT);
     673       31941 :   for (i=1; i<l; i++) gel(z,i) = QXQC_to_mod_shallow(gel(V,i), T);
     674        7476 :   return z;
     675             : }
     676             : 
     677             : GEN
     678     3826606 : RgX_renormalize_lg(GEN x, long lx)
     679             : {
     680             :   long i;
     681     6383140 :   for (i = lx-1; i>1; i--)
     682     5997491 :     if (! gequal0(gel(x,i))) break; /* _not_ isexactzero */
     683     3826606 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
     684     3826606 :   setlg(x, i+1); setsigne(x, i != 1); return x;
     685             : }
     686             : 
     687             : GEN
     688      534429 : RgV_to_RgX(GEN x, long v)
     689             : {
     690      534429 :   long i, k = lg(x);
     691             :   GEN p;
     692             : 
     693     1925086 :   while (--k && gequal0(gel(x,k)));
     694      534429 :   if (!k) return pol_0(v);
     695      527884 :   i = k+2; p = cgetg(i,t_POL);
     696      527884 :   p[1] = evalsigne(1) | evalvarn(v);
     697     5633190 :   x--; for (k=2; k<i; k++) gel(p,k) = gel(x,k);
     698      527884 :   return p;
     699             : }
     700             : GEN
     701      159428 : RgV_to_RgX_reverse(GEN x, long v)
     702             : {
     703      159428 :   long j, k, l = lg(x);
     704             :   GEN p;
     705             : 
     706      159428 :   for (k = 1; k < l; k++)
     707      159428 :     if (!gequal0(gel(x,k))) break;
     708      159428 :   if (k == l) return pol_0(v);
     709      159428 :   k -= 1;
     710      159428 :   l -= k;
     711      159428 :   x += k;
     712      159428 :   p = cgetg(l+1,t_POL);
     713      159428 :   p[1] = evalsigne(1) | evalvarn(v);
     714      910437 :   for (j=2, k=l; j<=l; j++) gel(p,j) = gel(x,--k);
     715      159428 :   return p;
     716             : }
     717             : 
     718             : /* return the (N-dimensional) vector of coeffs of p */
     719             : GEN
     720     8671773 : RgX_to_RgC(GEN x, long N)
     721             : {
     722             :   long i, l;
     723             :   GEN z;
     724     8671773 :   l = lg(x)-1; x++;
     725     8671773 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
     726     8671773 :   z = cgetg(N+1,t_COL);
     727    58017486 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
     728    16755384 :   for (   ; i<=N; i++) gel(z,i) = gen_0;
     729     8671853 :   return z;
     730             : }
     731             : GEN
     732      747000 : Rg_to_RgC(GEN x, long N)
     733             : {
     734      747000 :   return (typ(x) == t_POL)? RgX_to_RgC(x,N): scalarcol_shallow(x, N);
     735             : }
     736             : 
     737             : /* vector of polynomials (in v) whose coeffs are given by the columns of x */
     738             : GEN
     739       62050 : RgM_to_RgXV(GEN x, long v)
     740      346122 : { pari_APPLY_type(t_VEC, RgV_to_RgX(gel(x,i), v)) }
     741             : 
     742             : /* matrix whose entries are given by the coeffs of the polynomials in
     743             :  * vector v (considered as degree n-1 polynomials) */
     744             : GEN
     745      131748 : RgV_to_RgM(GEN x, long n)
     746      874948 : { pari_APPLY_type(t_MAT, Rg_to_RgC(gel(x,i), n)) }
     747             : 
     748             : GEN
     749        5797 : RgXV_to_RgM(GEN x, long n)
     750       36285 : { pari_APPLY_type(t_MAT, RgX_to_RgC(gel(x,i), n)) }
     751             : 
     752             : /* polynomial (in v) of polynomials (in w) whose coeffs are given by the columns of x */
     753             : GEN
     754       17656 : RgM_to_RgXX(GEN x, long v,long w)
     755             : {
     756       17656 :   long j, lx = lg(x);
     757       17656 :   GEN y = cgetg(lx+1, t_POL);
     758       17656 :   y[1] = evalsigne(1) | evalvarn(v);
     759       17656 :   y++;
     760      100560 :   for (j=1; j<lx; j++) gel(y,j) = RgV_to_RgX(gel(x,j), w);
     761       17656 :   return normalizepol_lg(--y, lx+1);
     762             : }
     763             : 
     764             : /* matrix whose entries are given by the coeffs of the polynomial v in
     765             :  * two variables (considered as degree n-1 polynomials) */
     766             : GEN
     767         308 : RgXX_to_RgM(GEN v, long n)
     768             : {
     769         308 :   long j, N = lg(v)-1;
     770         308 :   GEN y = cgetg(N, t_MAT);
     771        1001 :   for (j=1; j<N; j++) gel(y,j) = Rg_to_RgC(gel(v,j+1), n);
     772         308 :   return y;
     773             : }
     774             : 
     775             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     776             : GEN
     777       21660 : RgXY_swapspec(GEN x, long n, long w, long nx)
     778             : {
     779       21660 :   long j, ly = n+3;
     780       21660 :   GEN y = cgetg(ly, t_POL);
     781       21660 :   y[1] = evalsigne(1);
     782      266735 :   for (j=2; j<ly; j++)
     783             :   {
     784             :     long k;
     785      245075 :     GEN a = cgetg(nx+2,t_POL);
     786      245075 :     a[1] = evalsigne(1) | evalvarn(w);
     787     1236738 :     for (k=0; k<nx; k++)
     788             :     {
     789      991663 :       GEN xk = gel(x,k);
     790      991663 :       if (typ(xk)==t_POL)
     791      912728 :         gel(a,k+2) = j<lg(xk)? gel(xk,j): gen_0;
     792             :       else
     793       78935 :         gel(a,k+2) = j==2 ? xk: gen_0;
     794             :     }
     795      245075 :     gel(y,j) = normalizepol_lg(a, nx+2);
     796             :   }
     797       21660 :   return normalizepol_lg(y,ly);
     798             : }
     799             : 
     800             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     801             : GEN
     802         952 : RgXY_swap(GEN x, long n, long w)
     803             : {
     804         952 :   GEN z = RgXY_swapspec(x+2, n, w, lgpol(x));
     805         952 :   setvarn(z, varn(x)); return z;
     806             : }
     807             : 
     808             : long
     809         190 : RgXY_degreex(GEN b)
     810             : {
     811         190 :   long deg = 0, i;
     812         190 :   if (!signe(b)) return -1;
     813        1032 :   for (i = 2; i < lg(b); ++i)
     814             :   {
     815         842 :     GEN bi = gel(b, i);
     816         842 :     if (typ(bi) == t_POL)
     817         799 :       deg = maxss(deg, degpol(bi));
     818             :   }
     819         190 :   return deg;
     820             : }
     821             : 
     822             : /* return (x % X^n). Shallow */
     823             : GEN
     824     4562138 : RgXn_red_shallow(GEN a, long n)
     825             : {
     826     4562138 :   long i, L = n+2, l = lg(a);
     827             :   GEN  b;
     828     4562138 :   if (L >= l) return a; /* deg(x) < n */
     829     2971483 :   b = cgetg(L, t_POL); b[1] = a[1];
     830    23205157 :   for (i=2; i<L; i++) gel(b,i) = gel(a,i);
     831     2971483 :   return normalizepol_lg(b,L);
     832             : }
     833             : 
     834             : GEN
     835         357 : RgXnV_red_shallow(GEN x, long n)
     836        1680 : { pari_APPLY_type(t_VEC, RgXn_red_shallow(gel(x,i), n)) }
     837             : 
     838             : /* return (x * X^n). Shallow */
     839             : GEN
     840    74414832 : RgX_shift_shallow(GEN a, long n)
     841             : {
     842    74414832 :   long i, l = lg(a);
     843             :   GEN  b;
     844    74414832 :   if (l == 2 || !n) return a;
     845    54607709 :   l += n;
     846    54607709 :   if (n < 0)
     847             :   {
     848    44704023 :     if (l <= 2) return pol_0(varn(a));
     849    43915923 :     b = cgetg(l, t_POL); b[1] = a[1];
     850    43915923 :     a -= n;
     851   117520747 :     for (i=2; i<l; i++) gel(b,i) = gel(a,i);
     852             :   } else {
     853     9903686 :     b = cgetg(l, t_POL); b[1] = a[1];
     854     9903686 :     a -= n; n += 2;
     855    22665279 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     856   280614992 :     for (   ; i<l; i++) gel(b,i) = gel(a,i);
     857             :   }
     858    53819609 :   return b;
     859             : }
     860             : /* return (x * X^n). */
     861             : GEN
     862     1534405 : RgX_shift(GEN a, long n)
     863             : {
     864     1534405 :   long i, l = lg(a);
     865             :   GEN  b;
     866     1534405 :   if (l == 2 || !n) return RgX_copy(a);
     867     1534125 :   l += n;
     868     1534125 :   if (n < 0)
     869             :   {
     870         630 :     if (l <= 2) return pol_0(varn(a));
     871         588 :     b = cgetg(l, t_POL); b[1] = a[1];
     872         588 :     a -= n;
     873        2275 :     for (i=2; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     874             :   } else {
     875     1533495 :     b = cgetg(l, t_POL); b[1] = a[1];
     876     1533495 :     a -= n; n += 2;
     877     3859201 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     878     4277413 :     for (   ; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     879             :   }
     880     1534083 :   return b;
     881             : }
     882             : 
     883             : GEN
     884      316113 : RgX_rotate_shallow(GEN P, long k, long p)
     885             : {
     886      316113 :   long i, l = lgpol(P);
     887             :   GEN r;
     888      316113 :   if (signe(P)==0)
     889        1365 :     return pol_0(varn(P));
     890      314748 :   r = cgetg(p+2,t_POL); r[1] = P[1];
     891     2090452 :   for(i=0; i<p; i++)
     892             :   {
     893     1775704 :     long s = 2+(i+k)%p;
     894     1775704 :     gel(r,s) = i<l? gel(P,2+i): gen_0;
     895             :   }
     896      314748 :   return RgX_renormalize(r);
     897             : }
     898             : 
     899             : GEN
     900     2975178 : RgX_mulXn(GEN x, long d)
     901             : {
     902             :   pari_sp av;
     903             :   GEN z;
     904             :   long v;
     905     2975178 :   if (d >= 0) return RgX_shift(x, d);
     906     1463667 :   d = -d;
     907     1463667 :   v = RgX_val(x);
     908     1463667 :   if (v >= d) return RgX_shift(x, -d);
     909     1463653 :   av = avma;
     910     1463653 :   z = gred_rfrac_simple(RgX_shift_shallow(x, -v), pol_xn(d - v, varn(x)));
     911     1463653 :   return gerepileupto(av, z);
     912             : }
     913             : 
     914             : long
     915     2672610 : RgX_val(GEN x)
     916             : {
     917     2672610 :   long i, lx = lg(x);
     918     2672610 :   if (lx == 2) return LONG_MAX;
     919     2699455 :   for (i = 2; i < lx; i++)
     920     2699413 :     if (!isexactzero(gel(x,i))) break;
     921     2672470 :   if (i == lx) return LONG_MAX;/* possible with non-rational zeros */
     922     2672428 :   return i - 2;
     923             : }
     924             : long
     925    39107604 : RgX_valrem(GEN x, GEN *Z)
     926             : {
     927    39107604 :   long v, i, lx = lg(x);
     928    39107604 :   if (lx == 2) { *Z = pol_0(varn(x)); return LONG_MAX; }
     929    73939688 :   for (i = 2; i < lx; i++)
     930    73939689 :     if (!isexactzero(gel(x,i))) break;
     931             :   /* possible with non-rational zeros */
     932    39107606 :   if (i == lx) { *Z = pol_0(varn(x)); return LONG_MAX; }
     933    39107606 :   v = i - 2;
     934    39107606 :   *Z = RgX_shift_shallow(x, -v);
     935    39107607 :   return v;
     936             : }
     937             : long
     938       17227 : RgX_valrem_inexact(GEN x, GEN *Z)
     939             : {
     940             :   long v;
     941       17227 :   if (!signe(x)) { if (Z) *Z = pol_0(varn(x)); return LONG_MAX; }
     942       17220 :   for (v = 0;; v++)
     943       18389 :     if (!gequal0(gel(x,2+v))) break;
     944       17220 :   if (Z) *Z = RgX_shift_shallow(x, -v);
     945       17220 :   return v;
     946             : }
     947             : 
     948             : GEN
     949      109438 : RgXQC_red(GEN x, GEN T)
     950     2051574 : { pari_APPLY_type(t_COL, grem(gel(x,i), T)) }
     951             : 
     952             : GEN
     953        1148 : RgXQV_red(GEN x, GEN T)
     954       27545 : { pari_APPLY_type(t_VEC, grem(gel(x,i), T)) }
     955             : 
     956             : GEN
     957       14826 : RgXQM_red(GEN x, GEN T)
     958      124264 : { pari_APPLY_same(RgXQC_red(gel(x,i), T)) }
     959             : 
     960             : GEN
     961           0 : RgXQM_mul(GEN P, GEN Q, GEN T)
     962             : {
     963           0 :   return RgXQM_red(RgM_mul(P, Q), T);
     964             : }
     965             : 
     966             : GEN
     967      581458 : RgXQX_red(GEN P, GEN T)
     968             : {
     969      581458 :   long i, l = lg(P);
     970      581458 :   GEN Q = cgetg(l, t_POL);
     971      581458 :   Q[1] = P[1];
     972     3466453 :   for (i=2; i<l; i++) gel(Q,i) = grem(gel(P,i), T);
     973      581458 :   return normalizepol_lg(Q, l);
     974             : }
     975             : 
     976             : GEN
     977      231704 : RgX_deriv(GEN x)
     978             : {
     979      231704 :   long i,lx = lg(x)-1;
     980             :   GEN y;
     981             : 
     982      231704 :   if (lx<3) return pol_0(varn(x));
     983      230234 :   y = cgetg(lx,t_POL); gel(y,2) = gcopy(gel(x,3));
     984     1046795 :   for (i=3; i<lx ; i++) gel(y,i) = gmulsg(i-1,gel(x,i+1));
     985      230234 :   y[1] = x[1]; return normalizepol_lg(y,i);
     986             : }
     987             : 
     988             : GEN
     989      495559 : RgX_recipspec_shallow(GEN x, long l, long n)
     990             : {
     991             :   long i;
     992      495559 :   GEN z = cgetg(n+2,t_POL);
     993      495559 :   z[1] = 0; z += 2;
     994    15439534 :   for(i=0; i<l; i++)
     995    14943975 :     gel(z,n-i-1) = gel(x,i);
     996      623112 :   for(   ; i<n; i++)
     997      127553 :     gel(z, n-i-1) = gen_0;
     998      495559 :   return normalizepol_lg(z-2,n+2);
     999             : }
    1000             : 
    1001             : GEN
    1002       47732 : RgXn_recip_shallow(GEN P, long n)
    1003             : {
    1004       47732 :   GEN Q = RgX_recipspec_shallow(P+2, lgpol(P), n);
    1005       47732 :   setvarn(Q, varn(P));
    1006       47732 :   return Q;
    1007             : }
    1008             : 
    1009             : /* return coefficients s.t x = x_0 X^n + ... + x_n */
    1010             : GEN
    1011        9513 : RgX_recip(GEN x)
    1012             : {
    1013             :   long lx, i, j;
    1014        9513 :   GEN y = cgetg_copy(x, &lx);
    1015       92071 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gcopy(gel(x,j));
    1016        9513 :   return normalizepol_lg(y,lx);
    1017             : }
    1018             : /* shallow version */
    1019             : GEN
    1020      523744 : RgX_recip_shallow(GEN x)
    1021             : {
    1022             :   long lx, i, j;
    1023      523744 :   GEN y = cgetg_copy(x, &lx);
    1024     3618465 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
    1025      523752 :   return y;
    1026             : }
    1027             : /*******************************************************************/
    1028             : /*                                                                 */
    1029             : /*                      ADDITION / SUBTRACTION                     */
    1030             : /*                                                                 */
    1031             : /*******************************************************************/
    1032             : /* same variable */
    1033             : GEN
    1034    34091714 : RgX_add(GEN x, GEN y)
    1035             : {
    1036    34091714 :   long i, lx = lg(x), ly = lg(y);
    1037             :   GEN z;
    1038    34091714 :   if (ly <= lx) {
    1039    29819036 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1040   139788857 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1041    62983593 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1042    29819039 :     z = normalizepol_lg(z, lx);
    1043             :   } else {
    1044     4272678 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1045    28765869 :     for (i=2; i < lx; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1046    15149587 :     for (   ; i < ly; i++) gel(z,i) = gcopy(gel(y,i));
    1047     4272680 :     z = normalizepol_lg(z, ly);
    1048             :   }
    1049    34091718 :   return z;
    1050             : }
    1051             : GEN
    1052    15467704 : RgX_sub(GEN x, GEN y)
    1053             : {
    1054    15467704 :   long i, lx = lg(x), ly = lg(y);
    1055             :   GEN z;
    1056    15467704 :   if (ly <= lx) {
    1057    10474973 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1058    67185804 :     for (i=2; i < ly; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1059    19344978 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1060    10474966 :     z = normalizepol_lg(z, lx);
    1061             :   } else {
    1062     4992731 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1063    21074541 :     for (i=2; i < lx; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1064    13336489 :     for (   ; i < ly; i++) gel(z,i) = gneg(gel(y,i));
    1065     4992730 :     z = normalizepol_lg(z, ly);
    1066             :   }
    1067    15467704 :   return z;
    1068             : }
    1069             : GEN
    1070     1434683 : RgX_neg(GEN x)
    1071             : {
    1072     1434683 :   long i, lx = lg(x);
    1073     1434683 :   GEN y = cgetg(lx, t_POL); y[1] = x[1];
    1074    19862116 :   for (i=2; i<lx; i++) gel(y,i) = gneg(gel(x,i));
    1075     1434678 :   return y;
    1076             : }
    1077             : 
    1078             : GEN
    1079    12985667 : RgX_Rg_add(GEN y, GEN x)
    1080             : {
    1081             :   GEN z;
    1082    12985667 :   long lz = lg(y), i;
    1083    12985667 :   if (lz == 2) return scalarpol(x,varn(y));
    1084     9996411 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1085     9996411 :   gel(z,2) = gadd(gel(y,2),x);
    1086    36586432 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1087             :   /* probably useless unless lz = 3, but cannot be skipped if y is
    1088             :    * an inexact 0 */
    1089     9996411 :   return normalizepol_lg(z,lz);
    1090             : }
    1091             : GEN
    1092        2883 : RgX_Rg_add_shallow(GEN y, GEN x)
    1093             : {
    1094             :   GEN z;
    1095        2883 :   long lz = lg(y), i;
    1096        2883 :   if (lz == 2) return scalarpol(x,varn(y));
    1097        2883 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1098        2883 :   gel(z,2) = gadd(gel(y,2),x);
    1099        5836 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
    1100        2883 :   return z = normalizepol_lg(z,lz);
    1101             : }
    1102             : GEN
    1103       33179 : RgX_Rg_sub(GEN y, GEN x)
    1104             : {
    1105             :   GEN z;
    1106       33179 :   long lz = lg(y), i;
    1107       33179 :   if (lz == 2)
    1108             :   { /* scalarpol(gneg(x),varn(y)) optimized */
    1109        2576 :     long v = varn(y);
    1110        2576 :     if (isrationalzero(x)) return pol_0(v);
    1111          14 :     z = cgetg(3,t_POL);
    1112          14 :     z[1] = gequal0(x)? evalvarn(v)
    1113          14 :                    : evalvarn(v) | evalsigne(1);
    1114          14 :     gel(z,2) = gneg(x); return z;
    1115             :   }
    1116       30603 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1117       30603 :   gel(z,2) = gsub(gel(y,2),x);
    1118      131212 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1119       30603 :   return z = normalizepol_lg(z,lz);
    1120             : }
    1121             : GEN
    1122      597421 : Rg_RgX_sub(GEN x, GEN y)
    1123             : {
    1124             :   GEN z;
    1125      597421 :   long lz = lg(y), i;
    1126      597421 :   if (lz == 2) return scalarpol(x,varn(y));
    1127      596329 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1128      596329 :   gel(z,2) = gsub(x, gel(y,2));
    1129     1201291 :   for(i=3; i<lz; i++) gel(z,i) = gneg(gel(y,i));
    1130      596329 :   return z = normalizepol_lg(z,lz);
    1131             : }
    1132             : /*******************************************************************/
    1133             : /*                                                                 */
    1134             : /*                  KARATSUBA MULTIPLICATION                       */
    1135             : /*                                                                 */
    1136             : /*******************************************************************/
    1137             : #if 0
    1138             : /* to debug Karatsuba-like routines */
    1139             : GEN
    1140             : zx_debug_spec(GEN x, long nx)
    1141             : {
    1142             :   GEN z = cgetg(nx+2,t_POL);
    1143             :   long i;
    1144             :   for (i=0; i<nx; i++) gel(z,i+2) = stoi(x[i]);
    1145             :   z[1] = evalsigne(1); return z;
    1146             : }
    1147             : 
    1148             : GEN
    1149             : RgX_debug_spec(GEN x, long nx)
    1150             : {
    1151             :   GEN z = cgetg(nx+2,t_POL);
    1152             :   long i;
    1153             :   for (i=0; i<nx; i++) z[i+2] = x[i];
    1154             :   z[1] = evalsigne(1); return z;
    1155             : }
    1156             : #endif
    1157             : 
    1158             : /* generic multiplication */
    1159             : GEN
    1160     3408579 : RgX_addspec_shallow(GEN x, GEN y, long nx, long ny)
    1161             : {
    1162             :   GEN z, t;
    1163             :   long i;
    1164     3408579 :   if (nx == ny) {
    1165      729091 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1166     2740604 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1167      729091 :     return normalizepol_lg(z, nx+2);
    1168             :   }
    1169     2679488 :   if (ny < nx) {
    1170     2520177 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1171    11352529 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1172     7853734 :     for (   ; i < nx; i++) gel(t,i) = gel(x,i);
    1173     2520177 :     return normalizepol_lg(z, nx+2);
    1174             :   } else {
    1175      159311 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1176     3137647 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1177      400627 :     for (   ; i < ny; i++) gel(t,i) = gel(y,i);
    1178      159311 :     return normalizepol_lg(z, ny+2);
    1179             :   }
    1180             : }
    1181             : GEN
    1182      211212 : RgX_addspec(GEN x, GEN y, long nx, long ny)
    1183             : {
    1184             :   GEN z, t;
    1185             :   long i;
    1186      211212 :   if (nx == ny) {
    1187       12649 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1188     2143645 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1189       12649 :     return normalizepol_lg(z, nx+2);
    1190             :   }
    1191      198563 :   if (ny < nx) {
    1192      196862 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1193     3234843 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1194     2323379 :     for (   ; i < nx; i++) gel(t,i) = gcopy(gel(x,i));
    1195      196862 :     return normalizepol_lg(z, nx+2);
    1196             :   } else {
    1197        1701 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1198      327040 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1199       12040 :     for (   ; i < ny; i++) gel(t,i) = gcopy(gel(y,i));
    1200        1701 :     return normalizepol_lg(z, ny+2);
    1201             :   }
    1202             : }
    1203             : 
    1204             : /* Return the vector of coefficients of x, where we replace rational 0s by NULL
    1205             :  * [ to speed up basic operation s += x[i]*y[j] ]. We create a proper
    1206             :  * t_VECSMALL, to hold this, which can be left on stack: gerepile
    1207             :  * will not crash on it. The returned vector itself is not a proper GEN,
    1208             :  * we access the coefficients as x[i], i = 0..deg(x) */
    1209             : static GEN
    1210    28731137 : RgXspec_kill0(GEN x, long lx)
    1211             : {
    1212    28731137 :   GEN z = cgetg(lx+1, t_VECSMALL) + 1; /* inhibit gerepile-wise */
    1213             :   long i;
    1214   109281201 :   for (i=0; i <lx; i++)
    1215             :   {
    1216    80550064 :     GEN c = gel(x,i);
    1217    80550064 :     z[i] = (long)(isrationalzero(c)? NULL: c);
    1218             :   }
    1219    28731137 :   return z;
    1220             : }
    1221             : 
    1222             : INLINE GEN
    1223    56268918 : RgX_mulspec_basecase_limb(GEN x, GEN y, long a, long b)
    1224             : {
    1225    56268918 :   pari_sp av = avma;
    1226    56268918 :   GEN s = NULL;
    1227             :   long i;
    1228             : 
    1229   255477570 :   for (i=a; i<b; i++)
    1230   199208652 :     if (gel(y,i) && gel(x,-i))
    1231             :     {
    1232   119558988 :       GEN t = gmul(gel(y,i), gel(x,-i));
    1233   119558988 :       s = s? gadd(s, t): t;
    1234             :     }
    1235    56268918 :   return s? gerepileupto(av, s): gen_0;
    1236             : }
    1237             : 
    1238             : /* assume nx >= ny > 0, return x * y * t^v */
    1239             : static GEN
    1240    10218571 : RgX_mulspec_basecase(GEN x, GEN y, long nx, long ny, long v)
    1241             : {
    1242             :   long i, lz, nz;
    1243             :   GEN z;
    1244             : 
    1245    10218571 :   x = RgXspec_kill0(x,nx);
    1246    10218571 :   y = RgXspec_kill0(y,ny);
    1247    10218571 :   lz = nx + ny + 1; nz = lz-2;
    1248    10218571 :   lz += v;
    1249    10218571 :   z = cgetg(lz, t_POL) + 2; /* x:y:z [i] = term of degree i */
    1250    21975739 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1251    30480710 :   for (i=0; i<ny; i++)gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0, i+1);
    1252    19929714 :   for (  ; i<nx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ny);
    1253    20262139 :   for (  ; i<nz; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-nx+1,ny);
    1254    10218571 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1255             : }
    1256             : 
    1257             : /* return (x * X^d) + y. Assume d > 0 */
    1258             : GEN
    1259     2714090 : RgX_addmulXn_shallow(GEN x0, GEN y0, long d)
    1260             : {
    1261             :   GEN x, y, xd, yd, zd;
    1262             :   long a, lz, nx, ny;
    1263             : 
    1264     2714090 :   if (!signe(x0)) return y0;
    1265     2687189 :   ny = lgpol(y0);
    1266     2687189 :   nx = lgpol(x0);
    1267     2687189 :   zd = (GEN)avma;
    1268     2687189 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1269     2687189 :   if (a <= 0)
    1270             :   {
    1271      220064 :     lz = nx+d+2;
    1272      220064 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1273      876120 :     while (xd > x) gel(--zd,0) = gel(--xd,0);
    1274      220064 :     x = zd + a;
    1275      229995 :     while (zd > x) gel(--zd,0) = gen_0;
    1276             :   }
    1277             :   else
    1278             :   {
    1279     2467125 :     xd = new_chunk(d); yd = y+d;
    1280     2467125 :     x = RgX_addspec_shallow(x,yd, nx,a);
    1281     2467125 :     lz = (a>nx)? ny+2: lg(x)+d;
    1282    17186417 :     x += 2; while (xd > x) *--zd = *--xd;
    1283             :   }
    1284     8488708 :   while (yd > y) *--zd = *--yd;
    1285     2687189 :   *--zd = x0[1];
    1286     2687189 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1287             : }
    1288             : GEN
    1289      503484 : RgX_addmulXn(GEN x0, GEN y0, long d)
    1290             : {
    1291             :   GEN x, y, xd, yd, zd;
    1292             :   long a, lz, nx, ny;
    1293             : 
    1294      503484 :   if (!signe(x0)) return RgX_copy(y0);
    1295      502714 :   nx = lgpol(x0);
    1296      502714 :   ny = lgpol(y0);
    1297      502714 :   zd = (GEN)avma;
    1298      502714 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1299      502714 :   if (a <= 0)
    1300             :   {
    1301      291502 :     lz = nx+d+2;
    1302      291502 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1303     4259859 :     while (xd > x) gel(--zd,0) = gcopy(gel(--xd,0));
    1304      291502 :     x = zd + a;
    1305      757381 :     while (zd > x) gel(--zd,0) = gen_0;
    1306             :   }
    1307             :   else
    1308             :   {
    1309      211212 :     xd = new_chunk(d); yd = y+d;
    1310      211212 :     x = RgX_addspec(x,yd, nx,a);
    1311      211212 :     lz = (a>nx)? ny+2: lg(x)+d;
    1312     7842384 :     x += 2; while (xd > x) *--zd = *--xd;
    1313             :   }
    1314     2566390 :   while (yd > y) gel(--zd,0) = gcopy(gel(--yd,0));
    1315      502714 :   *--zd = x0[1];
    1316      502714 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1317             : }
    1318             : 
    1319             : /* return x * y mod t^n */
    1320             : static GEN
    1321     4143023 : RgXn_mul_basecase(GEN x, GEN y, long n)
    1322             : {
    1323     4143023 :   long i, lz = n+2, lx = lgpol(x), ly = lgpol(y);
    1324             :   GEN z;
    1325     4143023 :   if (lx < 0) return pol_0(varn(x));
    1326     4143023 :   if (ly < 0) return pol_0(varn(x));
    1327     4143023 :   z = cgetg(lz, t_POL) + 2;
    1328     4143023 :   x+=2; if (lx > n) lx = n;
    1329     4143023 :   y+=2; if (ly > n) ly = n;
    1330     4143023 :   z[-1] = x[-1];
    1331     4143023 :   if (ly > lx) { swap(x,y); lswap(lx,ly); }
    1332     4143023 :   x = RgXspec_kill0(x, lx);
    1333     4143023 :   y = RgXspec_kill0(y, ly);
    1334             :   /* x:y:z [i] = term of degree i */
    1335    18228515 :   for (i=0;i<ly; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,i+1);
    1336     6261720 :   for (  ; i<lx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ly);
    1337     4190902 :   for (  ; i<n; i++)  gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-lx+1,ly);
    1338     4143023 :   return normalizepol_lg(z - 2, lz);
    1339             : }
    1340             : /* Mulders / Karatsuba product f*g mod t^n (Hanrot-Zimmermann variant) */
    1341             : static GEN
    1342     5152273 : RgXn_mul2(GEN f, GEN g, long n)
    1343             : {
    1344     5152273 :   pari_sp av = avma;
    1345             :   GEN fe,fo, ge,go, l,h,m;
    1346             :   long n0, n1;
    1347     5152273 :   if (degpol(f) + degpol(g) < n) return RgX_mul(f,g);
    1348     4174425 :   if (n < 80) return RgXn_mul_basecase(f,g,n);
    1349       31402 :   n0 = n>>1; n1 = n-n0;
    1350       31402 :   RgX_even_odd(f, &fe, &fo);
    1351       31402 :   RgX_even_odd(g, &ge, &go);
    1352       31402 :   l = RgXn_mul(fe,ge,n1);
    1353       31402 :   h = RgXn_mul(fo,go,n0);
    1354       31402 :   m = RgX_sub(RgXn_mul(RgX_add(fe,fo),RgX_add(ge,go),n0), RgX_add(l,h));
    1355             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1356             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1357       31402 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1358             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1359       31402 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1360       31402 :   m = RgX_inflate(m,2);
    1361             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1362       31402 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1363       31402 :   h = RgX_inflate(h,2);
    1364       31402 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1365       31402 :   return gerepileupto(av, h);
    1366             : }
    1367             : /* (f*g) \/ x^n */
    1368             : static GEN
    1369      885463 : RgX_mulhigh_i2(GEN f, GEN g, long n)
    1370             : {
    1371      885463 :   long d = degpol(f)+degpol(g) + 1 - n;
    1372             :   GEN h;
    1373      885463 :   if (d <= 2) return RgX_shift_shallow(RgX_mul(f,g), -n);
    1374       21310 :   h = RgX_recip_shallow(RgXn_mul(RgX_recip_shallow(f),
    1375             :                                  RgX_recip_shallow(g), d));
    1376       21310 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1377             : }
    1378             : 
    1379             : /* (f*g) \/ x^n */
    1380             : static GEN
    1381           0 : RgX_sqrhigh_i2(GEN f, long n)
    1382             : {
    1383           0 :   long d = 2*degpol(f)+ 1 - n;
    1384             :   GEN h;
    1385           0 :   if (d <= 2) return RgX_shift_shallow(RgX_sqr(f), -n);
    1386           0 :   h = RgX_recip_shallow(RgXn_sqr(RgX_recip_shallow(f), d));
    1387           0 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1388             : }
    1389             : 
    1390             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
    1391             :  * b+2 were sent instead. na, nb = number of terms of a, b.
    1392             :  * Only c, c0, c1, c2 are genuine GEN.
    1393             :  */
    1394             : GEN
    1395    10930529 : RgX_mulspec(GEN a, GEN b, long na, long nb)
    1396             : {
    1397             :   GEN a0, c, c0;
    1398    10930529 :   long n0, n0a, i, v = 0;
    1399             :   pari_sp av;
    1400             : 
    1401    17332995 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v++; }
    1402    16819952 :   while (nb && isrationalzero(gel(b,0))) { b++; nb--; v++; }
    1403    10930529 :   if (na < nb) swapspec(a,b, na,nb);
    1404    10930529 :   if (!nb) return pol_0(0);
    1405             : 
    1406    10690614 :   if (nb < RgX_MUL_LIMIT) return RgX_mulspec_basecase(a,b,na,nb, v);
    1407      472043 :   RgX_shift_inplace_init(v);
    1408      472043 :   i = (na>>1); n0 = na-i; na = i;
    1409      472043 :   av = avma; a0 = a+n0; n0a = n0;
    1410     1327543 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1411             : 
    1412      472043 :   if (nb > n0)
    1413             :   {
    1414             :     GEN b0,c1,c2;
    1415             :     long n0b;
    1416             : 
    1417      470727 :     nb -= n0; b0 = b+n0; n0b = n0;
    1418     1438690 :     while (n0b && isrationalzero(gel(b,n0b-1))) n0b--;
    1419      470727 :     c = RgX_mulspec(a,b,n0a,n0b);
    1420      470727 :     c0 = RgX_mulspec(a0,b0, na,nb);
    1421             : 
    1422      470727 :     c2 = RgX_addspec_shallow(a0,a, na,n0a);
    1423      470727 :     c1 = RgX_addspec_shallow(b0,b, nb,n0b);
    1424             : 
    1425      470727 :     c1 = RgX_mulspec(c1+2,c2+2, lgpol(c1),lgpol(c2));
    1426      470727 :     c2 = RgX_sub(c1, RgX_add(c0,c));
    1427      470727 :     c0 = RgX_addmulXn_shallow(c0, c2, n0);
    1428             :   }
    1429             :   else
    1430             :   {
    1431        1316 :     c = RgX_mulspec(a,b,n0a,nb);
    1432        1316 :     c0 = RgX_mulspec(a0,b,na,nb);
    1433             :   }
    1434      472043 :   c0 = RgX_addmulXn(c0,c,n0);
    1435      472043 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1436             : }
    1437             : 
    1438             : INLINE GEN
    1439       41975 : RgX_sqrspec_basecase_limb(GEN x, long a, long i)
    1440             : {
    1441       41975 :   pari_sp av = avma;
    1442       41975 :   GEN s = NULL;
    1443       41975 :   long j, l = (i+1)>>1;
    1444      118174 :   for (j=a; j<l; j++)
    1445             :   {
    1446       76199 :     GEN xj = gel(x,j), xx = gel(x,i-j);
    1447       76199 :     if (xj && xx)
    1448             :     {
    1449       68177 :       GEN t = gmul(xj, xx);
    1450       68177 :       s = s? gadd(s, t): t;
    1451             :     }
    1452             :   }
    1453       41975 :   if (s) s = gshift(s,1);
    1454       41975 :   if ((i&1) == 0)
    1455             :   {
    1456       24962 :     GEN t = gel(x, i>>1);
    1457       24962 :     if (t) {
    1458       22365 :       t = gsqr(t);
    1459       22365 :       s = s? gadd(s, t): t;
    1460             :     }
    1461             :   }
    1462       41975 :   return s? gerepileupto(av,s): gen_0;
    1463             : }
    1464             : static GEN
    1465        7949 : RgX_sqrspec_basecase(GEN x, long nx, long v)
    1466             : {
    1467             :   long i, lz, nz;
    1468             :   GEN z;
    1469             : 
    1470        7949 :   if (!nx) return pol_0(0);
    1471        7949 :   x = RgXspec_kill0(x,nx);
    1472        7949 :   lz = (nx << 1) + 1, nz = lz-2;
    1473        7949 :   lz += v;
    1474        7949 :   z = cgetg(lz,t_POL) + 2;
    1475       17091 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1476       32911 :   for (i=0; i<nx; i++)gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1477       24962 :   for (  ; i<nz; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, i-nx+1, i);
    1478        7949 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1479             : }
    1480             : /* return x^2 mod t^n */
    1481             : static GEN
    1482           0 : RgXn_sqr_basecase(GEN x, long n)
    1483             : {
    1484           0 :   long i, lz = n+2, lx = lgpol(x);
    1485             :   GEN z;
    1486           0 :   if (lx < 0) return pol_0(varn(x));
    1487           0 :   z = cgetg(lz, t_POL);
    1488           0 :   z[1] = x[1];
    1489           0 :   x+=2; if (lx > n) lx = n;
    1490           0 :   x = RgXspec_kill0(x,lx);
    1491           0 :   z+=2;/* x:z [i] = term of degree i */
    1492           0 :   for (i=0;i<lx; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1493           0 :   for (  ; i<n; i++)  gel(z,i) = RgX_sqrspec_basecase_limb(x, i-lx+1, i);
    1494           0 :   z -= 2; return normalizepol_lg(z, lz);
    1495             : }
    1496             : /* Mulders / Karatsuba product f^2 mod t^n (Hanrot-Zimmermann variant) */
    1497             : static GEN
    1498         245 : RgXn_sqr2(GEN f, long n)
    1499             : {
    1500         245 :   pari_sp av = avma;
    1501             :   GEN fe,fo, l,h,m;
    1502             :   long n0, n1;
    1503         245 :   if (2*degpol(f) < n) return RgX_sqr_i(f);
    1504           0 :   if (n < 80) return RgXn_sqr_basecase(f,n);
    1505           0 :   n0 = n>>1; n1 = n-n0;
    1506           0 :   RgX_even_odd(f, &fe, &fo);
    1507           0 :   l = RgXn_sqr(fe,n1);
    1508           0 :   h = RgXn_sqr(fo,n0);
    1509           0 :   m = RgX_sub(RgXn_sqr(RgX_add(fe,fo),n0), RgX_add(l,h));
    1510             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1511             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1512           0 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1513             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1514           0 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1515           0 :   m = RgX_inflate(m,2);
    1516             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1517           0 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1518           0 :   h = RgX_inflate(h,2);
    1519           0 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1520           0 :   return gerepileupto(av, h);
    1521             : }
    1522             : GEN
    1523        7988 : RgX_sqrspec(GEN a, long na)
    1524             : {
    1525             :   GEN a0, c, c0, c1;
    1526        7988 :   long n0, n0a, i, v = 0;
    1527             :   pari_sp av;
    1528             : 
    1529       12559 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v += 2; }
    1530        7988 :   if (na<RgX_SQR_LIMIT) return RgX_sqrspec_basecase(a, na, v);
    1531          39 :   RgX_shift_inplace_init(v);
    1532          39 :   i = (na>>1); n0 = na-i; na = i;
    1533          39 :   av = avma; a0 = a+n0; n0a = n0;
    1534          39 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1535             : 
    1536          39 :   c = RgX_sqrspec(a,n0a);
    1537          39 :   c0 = RgX_sqrspec(a0,na);
    1538          39 :   c1 = gmul2n(RgX_mulspec(a0,a, na,n0a), 1);
    1539          39 :   c0 = RgX_addmulXn_shallow(c0,c1, n0);
    1540          39 :   c0 = RgX_addmulXn(c0,c,n0);
    1541          39 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1542             : }
    1543             : 
    1544             : /* (X^a + A)(X^b + B) - X^(a+b), where deg A < a, deg B < b */
    1545             : GEN
    1546      784587 : RgX_mul_normalized(GEN A, long a, GEN B, long b)
    1547             : {
    1548      784587 :   GEN z = RgX_mul(A, B);
    1549      784587 :   if (a < b)
    1550        5446 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(A, B, b-a), z, a);
    1551      779141 :   else if (a > b)
    1552      499149 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(B, A, a-b), z, b);
    1553             :   else
    1554      279992 :     z = RgX_addmulXn_shallow(RgX_add(A, B), z, a);
    1555      784587 :   return z;
    1556             : }
    1557             : 
    1558             : GEN
    1559     9515677 : RgX_mul_i(GEN x, GEN y)
    1560             : {
    1561     9515677 :   GEN z = RgX_mulspec(x+2, y+2, lgpol(x), lgpol(y));
    1562     9515677 :   setvarn(z, varn(x)); return z;
    1563             : }
    1564             : 
    1565             : GEN
    1566        7910 : RgX_sqr_i(GEN x)
    1567             : {
    1568        7910 :   GEN z = RgX_sqrspec(x+2, lgpol(x));
    1569        7910 :   setvarn(z,varn(x)); return z;
    1570             : }
    1571             : 
    1572             : /*******************************************************************/
    1573             : /*                                                                 */
    1574             : /*                               DIVISION                          */
    1575             : /*                                                                 */
    1576             : /*******************************************************************/
    1577             : GEN
    1578      577235 : RgX_Rg_divexact(GEN x, GEN y) {
    1579      577235 :   long i, lx = lg(x);
    1580             :   GEN z;
    1581      577235 :   if (lx == 2) return gcopy(x);
    1582      543089 :   switch(typ(y))
    1583             :   {
    1584      524133 :     case t_INT:
    1585      524133 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1586      490736 :       break;
    1587        3955 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1588             :   }
    1589      505737 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1590     4923110 :   for (i=2; i<lx; i++) gel(z,i) = gdivexact(gel(x,i),y);
    1591      505737 :   return z;
    1592             : }
    1593             : GEN
    1594    24677895 : RgX_Rg_div(GEN x, GEN y) {
    1595    24677895 :   long i, lx = lg(x);
    1596             :   GEN z;
    1597    24677895 :   if (lx == 2) return gcopy(x);
    1598    24459211 :   switch(typ(y))
    1599             :   {
    1600    15338933 :     case t_INT:
    1601    15338933 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1602     1123635 :       break;
    1603        4186 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1604             :   }
    1605    10239727 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1606    37776533 :   for (i=2; i<lx; i++) gel(z,i) = gdiv(gel(x,i),y);
    1607    10239724 :   return normalizepol_lg(z, lx);
    1608             : }
    1609             : GEN
    1610       13377 : RgX_normalize(GEN x)
    1611             : {
    1612       13377 :   GEN z, d = NULL;
    1613       13377 :   long i, n = lg(x)-1;
    1614       13377 :   for (i = n; i > 1; i--) { d = gel(x,i); if (!gequal0(d)) break; }
    1615       13377 :   if (i == 1) return pol_0(varn(x));
    1616       13377 :   if (i == n && isint1(d)) return x;
    1617        5019 :   n = i; z = cgetg(n+1, t_POL); z[1] = x[1];
    1618       11606 :   for (i=2; i<n; i++) gel(z,i) = gdiv(gel(x,i),d);
    1619        5019 :   gel(z,n) = Rg_get_1(d); return z;
    1620             : }
    1621             : GEN
    1622        2625 : RgX_divs(GEN x, long y) {
    1623             :   long i, lx;
    1624        2625 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1625       10318 :   for (i=2; i<lx; i++) gel(z,i) = gdivgs(gel(x,i),y);
    1626        2625 :   return normalizepol_lg(z, lx);
    1627             : }
    1628             : GEN
    1629       57550 : RgX_div_by_X_x(GEN a, GEN x, GEN *r)
    1630             : {
    1631       57550 :   long l = lg(a), i;
    1632             :   GEN a0, z0, z;
    1633             : 
    1634       57550 :   if (l <= 3)
    1635             :   {
    1636           0 :     if (r) *r = l == 2? gen_0: gcopy(gel(a,2));
    1637           0 :     return pol_0(0);
    1638             :   }
    1639       57550 :   z = cgetg(l-1, t_POL);
    1640       57550 :   z[1] = a[1];
    1641       57550 :   a0 = a + l-1;
    1642       57550 :   z0 = z + l-2; *z0 = *a0--;
    1643      907721 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1644             :   {
    1645      850171 :     GEN t = gadd(gel(a0--,0), gmul(x, gel(z0--,0)));
    1646      850171 :     gel(z0,0) = t;
    1647             :   }
    1648       57550 :   if (r) *r = gadd(gel(a0,0), gmul(x, gel(z0,0)));
    1649       57550 :   return z;
    1650             : }
    1651             : /* Polynomial division x / y:
    1652             :  *   if pr = ONLY_REM return remainder, otherwise return quotient
    1653             :  *   if pr = ONLY_DIVIDES return quotient if division is exact, else NULL
    1654             :  *   if pr != NULL set *pr to remainder, as the last object on stack */
    1655             : /* assume, typ(x) = typ(y) = t_POL, same variable */
    1656             : static GEN
    1657    12760870 : RgX_divrem_i(GEN x, GEN y, GEN *pr)
    1658             : {
    1659             :   pari_sp avy, av, av1;
    1660             :   long dx,dy,dz,i,j,sx,lr;
    1661             :   GEN z,p1,p2,rem,y_lead,mod,p;
    1662             :   GEN (*f)(GEN,GEN);
    1663             : 
    1664    12760870 :   if (!signe(y)) pari_err_INV("RgX_divrem",y);
    1665             : 
    1666    12760870 :   dy = degpol(y);
    1667    12760870 :   y_lead = gel(y,dy+2);
    1668    12760870 :   if (gequal0(y_lead)) /* normalize denominator if leading term is 0 */
    1669             :   {
    1670           0 :     pari_warn(warner,"normalizing a polynomial with 0 leading term");
    1671           0 :     for (dy--; dy>=0; dy--)
    1672             :     {
    1673           0 :       y_lead = gel(y,dy+2);
    1674           0 :       if (!gequal0(y_lead)) break;
    1675             :     }
    1676             :   }
    1677    12760870 :   if (!dy) /* y is constant */
    1678             :   {
    1679        3051 :     if (pr == ONLY_REM) return pol_0(varn(x));
    1680        3051 :     z = RgX_Rg_div(x, y_lead);
    1681        3051 :     if (pr == ONLY_DIVIDES) return z;
    1682        1168 :     if (pr) *pr = pol_0(varn(x));
    1683        1168 :     return z;
    1684             :   }
    1685    12757819 :   dx = degpol(x);
    1686    12757819 :   if (dx < dy)
    1687             :   {
    1688      695009 :     if (pr == ONLY_REM) return RgX_copy(x);
    1689      371761 :     if (pr == ONLY_DIVIDES) return signe(x)? NULL: pol_0(varn(x));
    1690      371740 :     z = pol_0(varn(x));
    1691      371740 :     if (pr) *pr = RgX_copy(x);
    1692      371740 :     return z;
    1693             :   }
    1694             : 
    1695             :   /* x,y in R[X], y non constant */
    1696    12062810 :   av = avma;
    1697    12062810 :   p = NULL;
    1698    12062810 :   if (RgX_is_FpX(x, &p) && RgX_is_FpX(y, &p) && p)
    1699             :   {
    1700      225841 :     z = FpX_divrem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p, pr);
    1701      225841 :     if (!z) return gc_NULL(av);
    1702      225841 :     z = FpX_to_mod(z, p);
    1703      225841 :     if (!pr || pr == ONLY_REM || pr == ONLY_DIVIDES)
    1704      121681 :       return gerepileupto(av, z);
    1705      104160 :     *pr = FpX_to_mod(*pr, p);
    1706      104160 :     gerepileall(av, 2, pr, &z);
    1707      104160 :     return z;
    1708             :   }
    1709    11836969 :   switch(typ(y_lead))
    1710             :   {
    1711           0 :     case t_REAL:
    1712           0 :       y_lead = ginv(y_lead);
    1713           0 :       f = gmul; mod = NULL;
    1714           0 :       break;
    1715        1880 :     case t_INTMOD:
    1716        1880 :     case t_POLMOD: y_lead = ginv(y_lead);
    1717        1880 :       f = gmul; mod = gmodulo(gen_1, gel(y_lead,1));
    1718        1880 :       break;
    1719    11835089 :     default: if (gequal1(y_lead)) y_lead = NULL;
    1720    11835089 :       f = gdiv; mod = NULL;
    1721             :   }
    1722             : 
    1723    11836969 :   if (y_lead == NULL)
    1724    10053644 :     p2 = gel(x,dx+2);
    1725             :   else {
    1726             :     for(;;) {
    1727     1783325 :       p2 = f(gel(x,dx+2),y_lead);
    1728     1783325 :       p2 = simplify_shallow(p2);
    1729     1783325 :       if (!isexactzero(p2) || (--dx < 0)) break;
    1730             :     }
    1731     1783325 :     if (dx < dy) /* leading coeff of x was in fact zero */
    1732             :     {
    1733           0 :       if (pr == ONLY_DIVIDES) {
    1734           0 :         set_avma(av);
    1735           0 :         return (dx < 0)? pol_0(varn(x)) : NULL;
    1736             :       }
    1737           0 :       if (pr == ONLY_REM)
    1738             :       {
    1739           0 :         if (dx < 0)
    1740           0 :           return gerepilecopy(av, scalarpol(p2, varn(x)));
    1741             :         else
    1742             :         {
    1743             :           GEN t;
    1744           0 :           set_avma(av);
    1745           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1746           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1747           0 :           return t;
    1748             :         }
    1749             :       }
    1750           0 :       if (pr) /* cf ONLY_REM above */
    1751             :       {
    1752           0 :         if (dx < 0)
    1753             :         {
    1754           0 :           p2 = gclone(p2);
    1755           0 :           set_avma(av);
    1756           0 :           z = pol_0(varn(x));
    1757           0 :           x = scalarpol(p2, varn(x));
    1758           0 :           gunclone(p2);
    1759             :         }
    1760             :         else
    1761             :         {
    1762             :           GEN t;
    1763           0 :           set_avma(av);
    1764           0 :           z = pol_0(varn(x));
    1765           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1766           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1767           0 :           x = t;
    1768             :         }
    1769           0 :         *pr = x;
    1770             :       }
    1771             :       else
    1772             :       {
    1773           0 :         set_avma(av);
    1774           0 :         z = pol_0(varn(x));
    1775             :       }
    1776           0 :       return z;
    1777             :     }
    1778             :   }
    1779             :   /* dx >= dy */
    1780    11836969 :   avy = avma;
    1781    11836969 :   dz = dx-dy;
    1782    11836969 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1783    11836969 :   x += 2;
    1784    11836969 :   z += 2;
    1785    11836969 :   y += 2;
    1786    11836969 :   gel(z,dz) = gcopy(p2);
    1787             : 
    1788    30847181 :   for (i=dx-1; i>=dy; i--)
    1789             :   {
    1790    19010212 :     av1=avma; p1=gel(x,i);
    1791  1068892932 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1792    19010212 :     if (y_lead) p1 = simplify(f(p1,y_lead));
    1793             : 
    1794    19010212 :     if (isrationalzero(p1)) { set_avma(av1); p1 = gen_0; }
    1795             :     else
    1796     9460216 :       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
    1797    19010212 :     gel(z,i-dy) = p1;
    1798             :   }
    1799    11836969 :   if (!pr) return gerepileupto(av,z-2);
    1800             : 
    1801     5149777 :   rem = (GEN)avma; av1 = (pari_sp)new_chunk(dx+3);
    1802     5149777 :   for (sx=0; ; i--)
    1803             :   {
    1804      349335 :     p1 = gel(x,i);
    1805             :     /* we always enter this loop at least once */
    1806    15169449 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1807     5499112 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1808     5499112 :     if (!gequal0(p1)) { sx = 1; break; } /* remainder is non-zero */
    1809     3421585 :     if (!isexactzero(p1)) break;
    1810     3405362 :     if (!i) break;
    1811      349335 :     set_avma(av1);
    1812             :   }
    1813     5149777 :   if (pr == ONLY_DIVIDES)
    1814             :   {
    1815        1659 :     if (sx) return gc_NULL(av);
    1816        1652 :     set_avma((pari_sp)rem); return gerepileupto(av,z-2);
    1817             :   }
    1818     5148118 :   lr=i+3; rem -= lr;
    1819     5148118 :   if (avma==av1) { set_avma((pari_sp)rem); p1 = gcopy(p1); }
    1820     5104548 :   else p1 = gerepileupto((pari_sp)rem,p1);
    1821     5148118 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1822     5148118 :   rem[1] = z[-1];
    1823     5148118 :   rem += 2;
    1824     5148118 :   gel(rem,i) = p1;
    1825     6617961 :   for (i--; i>=0; i--)
    1826             :   {
    1827     1469843 :     av1=avma; p1 = gel(x,i);
    1828     4666933 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1829     1469843 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1830     1469843 :     gel(rem,i) = avma==av1? gcopy(p1):gerepileupto(av1,p1);
    1831             :   }
    1832     5148118 :   rem -= 2;
    1833     5148118 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1834     5148118 :   if (pr == ONLY_REM) return gerepileupto(av,rem);
    1835     4149242 :   z -= 2;
    1836             :   {
    1837     4149242 :     GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
    1838     4149242 :     gerepilemanysp(av,avy,gptr,2); *pr = rem; return z;
    1839             :   }
    1840             : }
    1841             : 
    1842             : GEN
    1843    11438746 : RgX_divrem(GEN x, GEN y, GEN *pr)
    1844             : {
    1845    11438746 :   if (pr == ONLY_REM) return RgX_rem(x, y);
    1846    11438746 :   return RgX_divrem_i(x, y, pr);
    1847             : }
    1848             : 
    1849             : /* x and y in (R[Y]/T)[X]  (lifted), T in R[Y]. y preferably monic */
    1850             : GEN
    1851      111454 : RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)
    1852             : {
    1853             :   long vx, dx, dy, dz, i, j, sx, lr;
    1854             :   pari_sp av0, av, tetpil;
    1855             :   GEN z,p1,rem,lead;
    1856             : 
    1857      111454 :   if (!signe(y)) pari_err_INV("RgXQX_divrem",y);
    1858      111454 :   vx = varn(x);
    1859      111454 :   dx = degpol(x);
    1860      111454 :   dy = degpol(y);
    1861      111454 :   if (dx < dy)
    1862             :   {
    1863       34055 :     if (pr)
    1864             :     {
    1865       34055 :       av0 = avma; x = RgXQX_red(x, T);
    1866       34055 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: gen_0; }
    1867       34055 :       if (pr == ONLY_REM) return x;
    1868           0 :       *pr = x;
    1869             :     }
    1870           0 :     return pol_0(vx);
    1871             :   }
    1872       77399 :   lead = leading_coeff(y);
    1873       77399 :   if (!dy) /* y is constant */
    1874             :   {
    1875         546 :     if (pr && pr != ONLY_DIVIDES)
    1876             :     {
    1877           0 :       if (pr == ONLY_REM) return pol_0(vx);
    1878           0 :       *pr = pol_0(vx);
    1879             :     }
    1880         546 :     if (gequal1(lead)) return RgX_copy(x);
    1881           0 :     av0 = avma; x = gmul(x, ginvmod(lead,T)); tetpil = avma;
    1882           0 :     return gerepile(av0,tetpil,RgXQX_red(x,T));
    1883             :   }
    1884       76853 :   av0 = avma; dz = dx-dy;
    1885       76853 :   lead = gequal1(lead)? NULL: gclone(ginvmod(lead,T));
    1886       76853 :   set_avma(av0);
    1887       76853 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1888       76853 :   x += 2; y += 2; z += 2;
    1889             : 
    1890       76853 :   p1 = gel(x,dx); av = avma;
    1891       76853 :   gel(z,dz) = lead? gerepileupto(av, grem(gmul(p1,lead), T)): gcopy(p1);
    1892      339003 :   for (i=dx-1; i>=dy; i--)
    1893             :   {
    1894      262150 :     av=avma; p1=gel(x,i);
    1895     1548539 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1896      262150 :     if (lead) p1 = gmul(grem(p1, T), lead);
    1897      262150 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil, grem(p1, T));
    1898             :   }
    1899       76853 :   if (!pr) { guncloneNULL(lead); return z-2; }
    1900             : 
    1901       75607 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    1902       75607 :   for (sx=0; ; i--)
    1903             :   {
    1904       33199 :     p1 = gel(x,i);
    1905      391653 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1906      108806 :     tetpil=avma; p1 = grem(p1, T); if (!gequal0(p1)) { sx = 1; break; }
    1907       41606 :     if (!i) break;
    1908       33199 :     set_avma(av);
    1909             :   }
    1910       75607 :   if (pr == ONLY_DIVIDES)
    1911             :   {
    1912        5446 :     guncloneNULL(lead);
    1913        5446 :     if (sx) return gc_NULL(av0);
    1914        5313 :     set_avma((pari_sp)rem); return z-2;
    1915             :   }
    1916       70161 :   lr=i+3; rem -= lr;
    1917       70161 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1918       70161 :   rem[1] = z[-1];
    1919       70161 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    1920       70161 :   rem += 2; gel(rem,i) = p1;
    1921      153392 :   for (i--; i>=0; i--)
    1922             :   {
    1923       83231 :     av=avma; p1 = gel(x,i);
    1924      290747 :     for (j=0; j<=i && j<=dz; j++)
    1925      207516 :       p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1926       83231 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, grem(p1, T));
    1927             :   }
    1928       70161 :   rem -= 2;
    1929       70161 :   guncloneNULL(lead);
    1930       70161 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1931       70161 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    1932         126 :   *pr = rem; return z-2;
    1933             : }
    1934             : 
    1935             : /*******************************************************************/
    1936             : /*                                                                 */
    1937             : /*                        PSEUDO-DIVISION                          */
    1938             : /*                                                                 */
    1939             : /*******************************************************************/
    1940             : INLINE GEN
    1941      534513 : rem(GEN c, GEN T)
    1942             : {
    1943      534513 :   if (T && typ(c) == t_POL && varn(c) == varn(T)) c = RgX_rem(c, T);
    1944      534513 :   return c;
    1945             : }
    1946             : 
    1947             : /* x, y, are ZYX, lc(y) is an integer, T is a ZY */
    1948             : int
    1949        8451 : ZXQX_dvd(GEN x, GEN y, GEN T)
    1950             : {
    1951             :   long dx, dy, dz, i, p, T_ismonic;
    1952        8451 :   pari_sp av = avma, av2;
    1953             :   GEN y_lead;
    1954             : 
    1955        8451 :   if (!signe(y)) pari_err_INV("ZXQX_dvd",y);
    1956        8451 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1957        8451 :   if (typ(y_lead) == t_POL) y_lead = gel(y_lead, 2); /* t_INT */
    1958             :   /* if monic, no point in using pseudo-division */
    1959        8451 :   if (gequal1(y_lead)) return signe(RgXQX_rem(x, y, T)) == 0;
    1960        5357 :   T_ismonic = gequal1(leading_coeff(T));
    1961        5357 :   dx = degpol(x);
    1962        5357 :   if (dx < dy) return !signe(x);
    1963        5357 :   (void)new_chunk(2);
    1964        5357 :   x = RgX_recip_shallow(x)+2;
    1965        5357 :   y = RgX_recip_shallow(y)+2;
    1966             :   /* pay attention to sparse divisors */
    1967       11210 :   for (i = 1; i <= dy; i++)
    1968        5853 :     if (!signe(gel(y,i))) gel(y,i) = NULL;
    1969        5357 :   dz = dx-dy; p = dz+1;
    1970        5357 :   av2 = avma;
    1971             :   for (;;)
    1972       33669 :   {
    1973       39026 :     GEN m, x0 = gel(x,0), y0 = y_lead, cx = content(x0);
    1974       39026 :     x0 = gneg(x0); p--;
    1975       39026 :     m = gcdii(cx, y0);
    1976       39026 :     if (!equali1(m))
    1977             :     {
    1978       37701 :       x0 = gdiv(x0, m);
    1979       37701 :       y0 = diviiexact(y0, m);
    1980       37701 :       if (equali1(y0)) y0 = NULL;
    1981             :     }
    1982       81765 :     for (i=1; i<=dy; i++)
    1983             :     {
    1984       42739 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1985       42739 :       if (gel(y,i)) c = gadd(c, gmul(x0,gel(y,i)));
    1986       42739 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1987       42739 :       gel(x,i) = c;
    1988             :     }
    1989      415246 :     for (   ; i<=dx; i++)
    1990             :     {
    1991      376220 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1992      376220 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1993      376220 :       gel(x,i) = c;
    1994             :     }
    1995       44745 :     do { x++; dx--; } while (dx >= 0 && !signe(gel(x,0)));
    1996       39026 :     if (dx < dy) break;
    1997       33669 :     if (gc_needed(av2,1))
    1998             :     {
    1999           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZXQX_dvd dx = %ld >= %ld",dx,dy);
    2000           0 :       gerepilecoeffs(av2,x,dx+1);
    2001             :     }
    2002             :   }
    2003        5357 :   return gc_bool(av, dx < 0);
    2004             : }
    2005             : 
    2006             : /* T either NULL or a t_POL. */
    2007             : GEN
    2008       28763 : RgXQX_pseudorem(GEN x, GEN y, GEN T)
    2009             : {
    2010       28763 :   long vx = varn(x), dx, dy, dz, i, lx, p;
    2011       28763 :   pari_sp av = avma, av2;
    2012             :   GEN y_lead;
    2013             : 
    2014       28763 :   if (!signe(y)) pari_err_INV("RgXQX_pseudorem",y);
    2015       28763 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2016             :   /* if monic, no point in using pseudo-division */
    2017       28763 :   if (gequal1(y_lead)) return T? RgXQX_rem(x, y, T): RgX_rem(x, y);
    2018       25641 :   dx = degpol(x);
    2019       25641 :   if (dx < dy) return RgX_copy(x);
    2020       25634 :   (void)new_chunk(2);
    2021       25634 :   x = RgX_recip_shallow(x)+2;
    2022       25634 :   y = RgX_recip_shallow(y)+2;
    2023             :   /* pay attention to sparse divisors */
    2024       80458 :   for (i = 1; i <= dy; i++)
    2025       54824 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2026       25634 :   dz = dx-dy; p = dz+1;
    2027       25634 :   av2 = avma;
    2028             :   for (;;)
    2029             :   {
    2030       92484 :     gel(x,0) = gneg(gel(x,0)); p--;
    2031      299904 :     for (i=1; i<=dy; i++)
    2032             :     {
    2033      207420 :       GEN c = gmul(y_lead, gel(x,i));
    2034      207420 :       if (gel(y,i)) c = gadd(c, gmul(gel(x,0),gel(y,i)));
    2035      207420 :       gel(x,i) = rem(c, T);
    2036             :     }
    2037      319724 :     for (   ; i<=dx; i++)
    2038             :     {
    2039      227240 :       GEN c = gmul(y_lead, gel(x,i));
    2040      227240 :       gel(x,i) = rem(c, T);
    2041             :     }
    2042      102089 :     do { x++; dx--; } while (dx >= 0 && gequal0(gel(x,0)));
    2043       92484 :     if (dx < dy) break;
    2044       66850 :     if (gc_needed(av2,1))
    2045             :     {
    2046           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudorem dx = %ld >= %ld",dx,dy);
    2047           0 :       gerepilecoeffs(av2,x,dx+1);
    2048             :     }
    2049             :   }
    2050       25634 :   if (dx < 0) return pol_0(vx);
    2051       24101 :   lx = dx+3; x -= 2;
    2052       24101 :   x[0] = evaltyp(t_POL) | evallg(lx);
    2053       24101 :   x[1] = evalsigne(1) | evalvarn(vx);
    2054       24101 :   x = RgX_recip_shallow(x);
    2055       24101 :   if (p)
    2056             :   { /* multiply by y[0]^p   [beware dummy vars from FpX_FpXY_resultant] */
    2057        2655 :     GEN t = y_lead;
    2058        2655 :     if (T && typ(t) == t_POL && varn(t) == varn(T))
    2059           0 :       t = RgXQ_powu(t, p, T);
    2060             :     else
    2061        2655 :       t = gpowgs(t, p);
    2062        7299 :     for (i=2; i<lx; i++)
    2063             :     {
    2064        4644 :       GEN c = gmul(gel(x,i), t);
    2065        4644 :       gel(x,i) = rem(c,T);
    2066             :     }
    2067        2655 :     if (!T) return gerepileupto(av, x);
    2068             :   }
    2069       21446 :   return gerepilecopy(av, x);
    2070             : }
    2071             : 
    2072             : GEN
    2073       28763 : RgX_pseudorem(GEN x, GEN y) { return RgXQX_pseudorem(x,y, NULL); }
    2074             : 
    2075             : /* Compute z,r s.t lc(y)^(dx-dy+1) x = z y + r */
    2076             : GEN
    2077        7829 : RgXQX_pseudodivrem(GEN x, GEN y, GEN T, GEN *ptr)
    2078             : {
    2079        7829 :   long vx = varn(x), dx, dy, dz, i, iz, lx, lz, p;
    2080        7829 :   pari_sp av = avma, av2;
    2081             :   GEN z, r, ypow, y_lead;
    2082             : 
    2083        7829 :   if (!signe(y)) pari_err_INV("RgXQX_pseudodivrem",y);
    2084        7829 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2085        7829 :   if (gequal1(y_lead)) return T? RgXQX_divrem(x,y, T, ptr): RgX_divrem(x,y, ptr);
    2086        7257 :   dx = degpol(x);
    2087        7257 :   if (dx < dy) { *ptr = RgX_copy(x); return pol_0(vx); }
    2088        7257 :   if (dx == dy)
    2089             :   {
    2090          70 :     GEN x_lead = gel(x,lg(x)-1);
    2091          70 :     x = RgX_renormalize_lg(leafcopy(x), lg(x)-1);
    2092          70 :     y = RgX_renormalize_lg(leafcopy(y), lg(y)-1);
    2093          70 :     r = RgX_sub(RgX_Rg_mul(x, y_lead), RgX_Rg_mul(y, x_lead));
    2094          70 :     *ptr = gerepileupto(av, r); return scalarpol(x_lead, vx);
    2095             :   }
    2096        7187 :   (void)new_chunk(2);
    2097        7187 :   x = RgX_recip_shallow(x)+2;
    2098        7187 :   y = RgX_recip_shallow(y)+2;
    2099             :   /* pay attention to sparse divisors */
    2100       36794 :   for (i = 1; i <= dy; i++)
    2101       29607 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2102        7187 :   dz = dx-dy; p = dz+1;
    2103        7187 :   lz = dz+3;
    2104        7187 :   z = cgetg(lz, t_POL);
    2105        7187 :   z[1] = evalsigne(1) | evalvarn(vx);
    2106       27546 :   for (i = 2; i < lz; i++) gel(z,i) = gen_0;
    2107        7187 :   ypow = new_chunk(dz+1);
    2108        7187 :   gel(ypow,0) = gen_1;
    2109        7187 :   gel(ypow,1) = y_lead;
    2110       13172 :   for (i=2; i<=dz; i++)
    2111             :   {
    2112        5985 :     GEN c = gmul(gel(ypow,i-1), y_lead);
    2113        5985 :     gel(ypow,i) = rem(c,T);
    2114             :   }
    2115        7187 :   av2 = avma;
    2116        7187 :   for (iz=2;;)
    2117             :   {
    2118        7810 :     p--;
    2119       14997 :     gel(z,iz++) = rem(gmul(gel(x,0), gel(ypow,p)), T);
    2120       65083 :     for (i=1; i<=dy; i++)
    2121             :     {
    2122       50086 :       GEN c = gmul(y_lead, gel(x,i));
    2123       50086 :       if (gel(y,i)) c = gsub(c, gmul(gel(x,0),gel(y,i)));
    2124       50086 :       gel(x,i) = rem(c, T);
    2125             :     }
    2126       39138 :     for (   ; i<=dx; i++)
    2127             :     {
    2128       24141 :       GEN c = gmul(y_lead, gel(x,i));
    2129       24141 :       gel(x,i) = rem(c,T);
    2130             :     }
    2131       14997 :     x++; dx--;
    2132       20359 :     while (dx >= dy && gequal0(gel(x,0))) { x++; dx--; iz++; }
    2133       14997 :     if (dx < dy) break;
    2134        7810 :     if (gc_needed(av2,1))
    2135             :     {
    2136           0 :       GEN X = x-2;
    2137           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudodivrem dx=%ld >= %ld",dx,dy);
    2138           0 :       X[0] = evaltyp(t_POL)|evallg(dx+3); X[1] = z[1]; /* hack */
    2139           0 :       gerepileall(av2,2, &X, &z); x = X+2;
    2140             :     }
    2141             :   }
    2142       12941 :   while (dx >= 0 && gequal0(gel(x,0))) { x++; dx--; }
    2143        7187 :   if (dx < 0)
    2144          98 :     x = pol_0(vx);
    2145             :   else
    2146             :   {
    2147        7089 :     lx = dx+3; x -= 2;
    2148        7089 :     x[0] = evaltyp(t_POL) | evallg(lx);
    2149        7089 :     x[1] = evalsigne(1) | evalvarn(vx);
    2150        7089 :     x = RgX_recip_shallow(x);
    2151             :   }
    2152        7187 :   z = RgX_recip_shallow(z);
    2153        7187 :   r = x;
    2154        7187 :   if (p)
    2155             :   {
    2156        3234 :     GEN c = gel(ypow,p); r = RgX_Rg_mul(r, c);
    2157        3234 :     if (T && typ(c) == t_POL && varn(c) == varn(T)) r = RgXQX_red(r, T);
    2158             :   }
    2159        7187 :   gerepileall(av, 2, &z, &r);
    2160        7187 :   *ptr = r; return z;
    2161             : }
    2162             : GEN
    2163        7626 : RgX_pseudodivrem(GEN x, GEN y, GEN *ptr)
    2164        7626 : { return RgXQX_pseudodivrem(x,y,NULL,ptr); }
    2165             : 
    2166             : GEN
    2167       12789 : RgXQX_mul(GEN x, GEN y, GEN T)
    2168             : {
    2169       12789 :   return RgXQX_red(RgX_mul(x,y), T);
    2170             : }
    2171             : GEN
    2172    74756392 : RgX_Rg_mul(GEN y, GEN x) {
    2173             :   long i, ly;
    2174    74756392 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2175    74756387 :   if (ly == 2) return z;
    2176   284903548 :   for (i = 2; i < ly; i++) gel(z,i) = gmul(x,gel(y,i));
    2177    74669005 :   return normalizepol_lg(z,ly);
    2178             : }
    2179             : GEN
    2180         539 : RgX_muls(GEN y, long x) {
    2181             :   long i, ly;
    2182         539 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2183         539 :   if (ly == 2) return z;
    2184        3192 :   for (i = 2; i < ly; i++) gel(z,i) = gmulsg(x,gel(y,i));
    2185         448 :   return normalizepol_lg(z,ly);
    2186             : }
    2187             : GEN
    2188        4417 : RgXQX_RgXQ_mul(GEN x, GEN y, GEN T)
    2189             : {
    2190        4417 :   return RgXQX_red(RgX_Rg_mul(x,y), T);
    2191             : }
    2192             : GEN
    2193          77 : RgXQV_RgXQ_mul(GEN v, GEN x, GEN T)
    2194             : {
    2195          77 :   return RgXQV_red(RgV_Rg_mul(v,x), T);
    2196             : }
    2197             : 
    2198             : GEN
    2199        1792 : RgXQX_sqr(GEN x, GEN T)
    2200             : {
    2201        1792 :   return RgXQX_red(RgX_sqr(x), T);
    2202             : }
    2203             : 
    2204             : GEN
    2205         448 : RgXQX_powers(GEN P, long n, GEN T)
    2206             : {
    2207         448 :   GEN v = cgetg(n+2, t_VEC);
    2208             :   long i;
    2209         448 :   gel(v, 1) = pol_1(varn(T));
    2210         448 :   if (n==0) return v;
    2211         448 :   gel(v, 2) = gcopy(P);
    2212        1512 :   for (i = 2; i <= n; i++) gel(v,i+1) = RgXQX_mul(P, gel(v,i), T);
    2213         448 :   return v;
    2214             : }
    2215             : 
    2216             : static GEN
    2217      141053 : _add(void *data, GEN x, GEN y) { (void)data; return RgX_add(x, y); }
    2218             : static GEN
    2219           0 : _sub(void *data, GEN x, GEN y) { (void)data; return RgX_sub(x, y); }
    2220             : static GEN
    2221       57303 : _sqr(void *data, GEN x) { return RgXQ_sqr(x, (GEN)data); }
    2222             : static GEN
    2223       47985 : _mul(void *data, GEN x, GEN y) { return RgXQ_mul(x,y, (GEN)data); }
    2224             : static GEN
    2225      253537 : _cmul(void *data, GEN P, long a, GEN x) { (void)data; return RgX_Rg_mul(x,gel(P,a+2)); }
    2226             : static GEN
    2227      251107 : _one(void *data) { return pol_1(varn((GEN)data)); }
    2228             : static GEN
    2229         462 : _zero(void *data) { return pol_0(varn((GEN)data)); }
    2230             : static GEN
    2231      152909 : _red(void *data, GEN x) { (void)data; return gcopy(x); }
    2232             : 
    2233             : static struct bb_algebra RgXQ_algebra = { _red, _add, _sub,
    2234             :               _mul, _sqr, _one, _zero };
    2235             : 
    2236             : GEN
    2237           0 : RgX_RgXQV_eval(GEN Q, GEN x, GEN T)
    2238             : {
    2239           0 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)T,&RgXQ_algebra,_cmul);
    2240             : }
    2241             : 
    2242             : GEN
    2243      112225 : RgX_RgXQ_eval(GEN Q, GEN x, GEN T)
    2244             : {
    2245      112225 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2246      112225 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)T,&RgXQ_algebra,_cmul);
    2247             : }
    2248             : 
    2249             : /* mod X^n */
    2250             : struct modXn {
    2251             :   long v; /* varn(X) */
    2252             :   long n;
    2253             : } ;
    2254             : static GEN
    2255       11613 : _sqrXn(void *data, GEN x) {
    2256       11613 :   struct modXn *S = (struct modXn*)data;
    2257       11613 :   return RgXn_sqr(x, S->n);
    2258             : }
    2259             : static GEN
    2260        4234 : _mulXn(void *data, GEN x, GEN y) {
    2261        4234 :   struct modXn *S = (struct modXn*)data;
    2262        4234 :   return RgXn_mul(x,y, S->n);
    2263             : }
    2264             : static GEN
    2265        1477 : _oneXn(void *data) {
    2266        1477 :   struct modXn *S = (struct modXn*)data;
    2267        1477 :   return pol_1(S->v);
    2268             : }
    2269             : static GEN
    2270           0 : _zeroXn(void *data) {
    2271           0 :   struct modXn *S = (struct modXn*)data;
    2272           0 :   return pol_0(S->v);
    2273             : }
    2274             : static struct bb_algebra RgXn_algebra = { _red, _add, _sub, _mulXn, _sqrXn,
    2275             :                                           _oneXn, _zeroXn };
    2276             : 
    2277             : GEN
    2278         357 : RgXn_powers(GEN x, long m, long n)
    2279             : {
    2280         357 :   long d = degpol(x);
    2281         357 :   int use_sqr = (d<<1) >= n;
    2282             :   struct modXn S;
    2283         357 :   S.v = varn(x); S.n = n;
    2284         357 :   return gen_powers(x,m,use_sqr,(void*)&S,_sqrXn,_mulXn,_oneXn);
    2285             : }
    2286             : 
    2287             : GEN
    2288        2216 : RgXn_powu_i(GEN x, ulong m, long n)
    2289             : {
    2290             :   struct modXn S;
    2291             :   long v;
    2292        2216 :   if (n == 1) return x;
    2293        2188 :   v = RgX_valrem(x, &x);
    2294        2188 :   if (v) n -= m * v;
    2295        2188 :   S.v = varn(x); S.n = n;
    2296        2188 :   x = gen_powu_i(x, m, (void*)&S,_sqrXn,_mulXn);
    2297        2188 :   if (v) x = RgX_shift_shallow(x, m * v);
    2298        2188 :   return x;
    2299             : }
    2300             : GEN
    2301           0 : RgXn_powu(GEN x, ulong m, long n)
    2302             : {
    2303             :   pari_sp av;
    2304           0 :   if (n == 1) return gcopy(x);
    2305           0 :   av = avma; return gerepilecopy(av, RgXn_powu_i(x, m, n));
    2306             : }
    2307             : 
    2308             : GEN
    2309         714 : RgX_RgXnV_eval(GEN Q, GEN x, long n)
    2310             : {
    2311             :   struct modXn S;
    2312         714 :   S.v = varn(gel(x,2)); S.n = n;
    2313         714 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)&S,&RgXn_algebra,_cmul);
    2314             : }
    2315             : 
    2316             : GEN
    2317           0 : RgX_RgXn_eval(GEN Q, GEN x, long n)
    2318             : {
    2319           0 :   int use_sqr = 2*degpol(x) >= n;
    2320             :   struct modXn S;
    2321           0 :   S.v = varn(x); S.n = n;
    2322           0 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2323             : }
    2324             : 
    2325             : /* Q(x) mod t^n, x in R[t], n >= 1 */
    2326             : GEN
    2327        3388 : RgXn_eval(GEN Q, GEN x, long n)
    2328             : {
    2329        3388 :   long d = degpol(x);
    2330             :   int use_sqr;
    2331             :   struct modXn S;
    2332        3388 :   if (d == 1 && isrationalzero(gel(x,2)))
    2333             :   {
    2334        3381 :     GEN y = RgX_unscale(Q, gel(x,3));
    2335        3381 :     setvarn(y, varn(x)); return y;
    2336             :   }
    2337           7 :   S.v = varn(x);
    2338           7 :   S.n = n;
    2339           7 :   use_sqr = (d<<1) >= n;
    2340           7 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2341             : }
    2342             : 
    2343             : /* (f*g mod t^n) \ t^n2, assuming 2*n2 >= n */
    2344             : static GEN
    2345     1475132 : RgXn_mulhigh(GEN f, GEN g, long n2, long n)
    2346             : {
    2347     1475132 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2348     1475132 :   return RgX_add(RgX_mulhigh_i(fl, g, n2), RgXn_mul(fh, g, n - n2));
    2349             : }
    2350             : 
    2351             : /* (f^2 mod t^n) \ t^n2, assuming 2*n2 >= n */
    2352             : static GEN
    2353           0 : RgXn_sqrhigh(GEN f, long n2, long n)
    2354             : {
    2355           0 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2356           0 :   return RgX_add(RgX_mulhigh_i(fl, f, n2), RgXn_mul(fh, f, n - n2));
    2357             : }
    2358             : 
    2359             : GEN
    2360      904806 : RgXn_inv_i(GEN f, long e)
    2361             : {
    2362             :   pari_sp av;
    2363             :   ulong mask;
    2364             :   GEN W, a;
    2365      904806 :   long v = varn(f), n = 1;
    2366             : 
    2367      904806 :   if (!signe(f)) pari_err_INV("RgXn_inv",f);
    2368      904799 :   a = ginv(gel(f,2));
    2369      904789 :   if (e == 1) return scalarpol(a, v);
    2370      901500 :   else if (e == 2)
    2371             :   {
    2372             :     GEN b;
    2373      214371 :     if (degpol(f) <= 0 || gequal0(b = gel(f,3))) return scalarpol(a, v);
    2374      192677 :     av = avma; b = gneg(b);
    2375      192680 :     if (!gequal1(a)) b = gmul(b, gsqr(a));
    2376      192682 :     W = deg1pol_shallow(b, a, v);
    2377      192683 :     return gcopy(W);
    2378             :   }
    2379      687129 :   W = scalarpol_shallow(a,v);
    2380      687129 :   mask = quadratic_prec_mask(e);
    2381      687129 :   av = avma;
    2382     2161022 :   for (;mask>1;)
    2383             :   {
    2384             :     GEN u, fr;
    2385     1473893 :     long n2 = n;
    2386     1473893 :     n<<=1; if (mask & 1) n--;
    2387     1473893 :     mask >>= 1;
    2388     1473893 :     fr = RgXn_red_shallow(f, n);
    2389     1473893 :     u = RgXn_mul(W, RgXn_mulhigh(fr, W, n2, n), n-n2);
    2390     1473893 :     W = RgX_sub(W, RgX_shift_shallow(u, n2));
    2391     1473893 :     if (gc_needed(av,2))
    2392             :     {
    2393           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_inv, e = %ld", n);
    2394           0 :       W = gerepileupto(av, W);
    2395             :     }
    2396             :   }
    2397      687129 :   return W;
    2398             : }
    2399             : 
    2400             : static GEN
    2401           0 : RgXn_inv_FpX(GEN x, long e, GEN p)
    2402             : {
    2403           0 :   pari_sp av = avma;
    2404             :   GEN r;
    2405           0 :   if (lgefint(p) == 3)
    2406             :   {
    2407           0 :     ulong pp = uel(p, 2);
    2408           0 :     if (pp == 2)
    2409           0 :       r = F2x_to_ZX(F2xn_inv(RgX_to_F2x(x), e));
    2410             :     else
    2411           0 :       r = Flx_to_ZX_inplace(Flxn_inv(RgX_to_Flx(x, pp), e, pp));
    2412             :   }
    2413             :   else
    2414           0 :     r = FpXn_inv(RgX_to_FpX(x, p), e, p);
    2415           0 :   return gerepileupto(av, FpX_to_mod(r, p));
    2416             : }
    2417             : 
    2418             : static GEN
    2419           0 : RgXn_inv_FpXQX(GEN x, long n, GEN pol, GEN p)
    2420             : {
    2421           0 :   pari_sp av = avma;
    2422           0 :   GEN r, T = RgX_to_FpX(pol, p);
    2423           0 :   if (signe(T) == 0) pari_err_OP("/", gen_1, x);
    2424           0 :   r = FpXQXn_inv(RgX_to_FpXQX(x, T, p), n, T, p);
    2425           0 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2426             : }
    2427             : 
    2428             : #define code(t1,t2) ((t1 << 6) | t2)
    2429             : 
    2430             : static GEN
    2431       77904 : RgXn_inv_fast(GEN x, long e)
    2432             : {
    2433             :   GEN p, pol;
    2434             :   long pa;
    2435       77904 :   long t = RgX_type(x,&p,&pol,&pa);
    2436       77903 :   switch(t)
    2437             :   {
    2438           0 :     case t_INTMOD: return RgXn_inv_FpX(x, e, p);
    2439           0 :     case code(t_POLMOD, t_INTMOD):
    2440           0 :                    return RgXn_inv_FpXQX(x, e, pol, p);
    2441       77903 :     default:       return NULL;
    2442             :   }
    2443             : }
    2444             : #undef code
    2445             : 
    2446             : GEN
    2447       77904 : RgXn_inv(GEN f, long e)
    2448             : {
    2449       77904 :   GEN h = RgXn_inv_fast(f, e);
    2450       77903 :   if (h) return h;
    2451       77903 :   return RgXn_inv_i(f, e);
    2452             : }
    2453             : 
    2454             : /* Compute intformal(x^n*S)/x^(n+1) */
    2455             : static GEN
    2456       14210 : RgX_integXn(GEN x, long n)
    2457             : {
    2458       14210 :   long i, lx = lg(x);
    2459             :   GEN y;
    2460       14210 :   if (lx == 2) return RgX_copy(x);
    2461       13846 :   y = cgetg(lx, t_POL); y[1] = x[1];
    2462       27713 :   for (i=2; i<lx; i++)
    2463       13867 :     gel(y,i) = gdivgs(gel(x,i), n+i-1);
    2464       13846 :   return RgX_renormalize_lg(y, lx);;
    2465             : }
    2466             : 
    2467             : GEN
    2468       12971 : RgXn_expint(GEN h, long e)
    2469             : {
    2470       12971 :   pari_sp av = avma, av2;
    2471       12971 :   long v = varn(h), n;
    2472       12971 :   GEN f = pol_1(v), g;
    2473             :   ulong mask;
    2474             : 
    2475       12971 :   if (!signe(h)) return f;
    2476       12971 :   g = pol_1(v);
    2477       12971 :   n = 1; mask = quadratic_prec_mask(e);
    2478       12971 :   av2 = avma;
    2479       14210 :   for (;mask>1;)
    2480             :   {
    2481             :     GEN u, w;
    2482       14210 :     long n2 = n;
    2483       14210 :     n<<=1; if (mask & 1) n--;
    2484       14210 :     mask >>= 1;
    2485       14210 :     u = RgXn_mul(g, RgX_mulhigh_i(f, RgXn_red_shallow(h, n2-1), n2-1), n-n2);
    2486       14210 :     u = RgX_add(u, RgX_shift_shallow(RgXn_red_shallow(h, n-1), 1-n2));
    2487       14210 :     w = RgXn_mul(f, RgX_integXn(u, n2-1), n-n2);
    2488       14210 :     f = RgX_add(f, RgX_shift_shallow(w, n2));
    2489       14210 :     if (mask<=1) break;
    2490        1239 :     u = RgXn_mul(g, RgXn_mulhigh(f, g, n2, n), n-n2);
    2491        1239 :     g = RgX_sub(g, RgX_shift_shallow(u, n2));
    2492        1239 :     if (gc_needed(av2,2))
    2493             :     {
    2494           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"RgXn_expint, e = %ld", n);
    2495           0 :       gerepileall(av2, 2, &f, &g);
    2496             :     }
    2497             :   }
    2498       12971 :   return gerepileupto(av, f);
    2499             : }
    2500             : 
    2501             : GEN
    2502           0 : RgXn_exp(GEN h, long e)
    2503             : {
    2504           0 :   long d = degpol(h);
    2505           0 :   if (d < 0) return pol_1(varn(h));
    2506           0 :   if (!d || !gequal0(gel(h,2)))
    2507           0 :     pari_err_DOMAIN("RgXn_exp","valuation", "<", gen_1, h);
    2508           0 :   return RgXn_expint(RgX_deriv(h), e);
    2509             : }
    2510             : 
    2511             : GEN
    2512          91 : RgXn_reverse(GEN f, long e)
    2513             : {
    2514          91 :   pari_sp av = avma, av2;
    2515             :   ulong mask;
    2516             :   GEN fi, a, df, W, an;
    2517          91 :   long v = varn(f), n=1;
    2518          91 :   if (degpol(f)<1 || !gequal0(gel(f,2)))
    2519           0 :     pari_err_INV("serreverse",f);
    2520          91 :   fi = ginv(gel(f,3));
    2521          91 :   a = deg1pol_shallow(fi,gen_0,v);
    2522          91 :   if (e <= 2) return gerepilecopy(av, a);
    2523          91 :   W = scalarpol(fi,v);
    2524          91 :   df = RgX_deriv(f);
    2525          91 :   mask = quadratic_prec_mask(e);
    2526          91 :   av2 = avma;
    2527         448 :   for (;mask>1;)
    2528             :   {
    2529             :     GEN u, fa, fr;
    2530         357 :     long n2 = n, rt;
    2531         357 :     n<<=1; if (mask & 1) n--;
    2532         357 :     mask >>= 1;
    2533         357 :     fr = RgXn_red_shallow(f, n);
    2534         357 :     rt = brent_kung_optpow(degpol(fr), 4, 3);
    2535         357 :     an = RgXn_powers(a, rt, n);
    2536         357 :     if (n>1)
    2537             :     {
    2538         357 :       long n4 = (n2+1)>>1;
    2539         357 :       GEN dfr = RgXn_red_shallow(df, n2);
    2540         357 :       dfr = RgX_RgXnV_eval(dfr, RgXnV_red_shallow(an, n2), n2);
    2541         357 :       u = RgX_shift(RgX_Rg_sub(RgXn_mul(W, dfr, n2), gen_1), -n4);
    2542         357 :       W = RgX_sub(W, RgX_shift(RgXn_mul(u, W, n2-n4), n4));
    2543             :     }
    2544         357 :     fa = RgX_sub(RgX_RgXnV_eval(fr, an, n), pol_x(v));
    2545         357 :     fa = RgX_shift(fa, -n2);
    2546         357 :     a = RgX_sub(a, RgX_shift(RgXn_mul(W, fa, n-n2), n2));
    2547         357 :     if (gc_needed(av2,2))
    2548             :     {
    2549           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_reverse, e = %ld", n);
    2550           0 :       gerepileall(av2, 2, &a, &W);
    2551             :     }
    2552             :   }
    2553          91 :   return gerepileupto(av, a);
    2554             : }
    2555             : 
    2556             : GEN
    2557           0 : RgXn_sqrt(GEN h, long e)
    2558             : {
    2559           0 :   pari_sp av = avma, av2;
    2560           0 :   long v = varn(h), n = 1;
    2561           0 :   GEN f = scalarpol(gen_1, v), df = f;
    2562           0 :   ulong mask = quadratic_prec_mask(e);
    2563           0 :   if (degpol(h)<0 || !gequal1(gel(h,2)))
    2564           0 :     pari_err_SQRTN("RgXn_sqrt",h);
    2565           0 :   av2 = avma;
    2566             :   while(1)
    2567           0 :   {
    2568           0 :     long n2 = n, m;
    2569             :     GEN g;
    2570           0 :     n<<=1; if (mask & 1) n--;
    2571           0 :     mask >>= 1;
    2572           0 :     m = n-n2;
    2573           0 :     g = RgX_sub(RgXn_sqrhigh(f, n2, n), RgX_shift_shallow(RgXn_red_shallow(h, n),-n2));
    2574           0 :     f = RgX_sub(f, RgX_shift_shallow(RgXn_mul(gmul2n(df, -1), g, m), n2));
    2575           0 :     if (mask==1) return gerepileupto(av, f);
    2576           0 :     g = RgXn_mul(df, RgXn_mulhigh(df, f, n2, n), m);
    2577           0 :     df = RgX_sub(df, RgX_shift_shallow(g, n2));
    2578           0 :     if (gc_needed(av2,2))
    2579             :     {
    2580           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_sqrt, e = %ld", n);
    2581           0 :       gerepileall(av2, 2, &f, &df);
    2582             :     }
    2583             :   }
    2584             : }
    2585             : 
    2586             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2587             : GEN
    2588       19063 : RgXQ_powu(GEN x, ulong n, GEN T)
    2589             : {
    2590       19063 :   pari_sp av = avma;
    2591             : 
    2592       19063 :   if (!n) return pol_1(varn(x));
    2593       17228 :   if (n == 1) return RgX_copy(x);
    2594       15281 :   x = gen_powu_i(x, n, (void*)T, &_sqr, &_mul);
    2595       15281 :   return gerepilecopy(av, x);
    2596             : }
    2597             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2598             : GEN
    2599       19108 : RgXQ_pow(GEN x, GEN n, GEN T)
    2600             : {
    2601             :   pari_sp av;
    2602       19108 :   long s = signe(n);
    2603             : 
    2604       19108 :   if (!s) return pol_1(varn(x));
    2605       19108 :   if (is_pm1(n) == 1)
    2606           0 :     return (s < 0)? RgXQ_inv(x, T): RgX_copy(x);
    2607       19108 :   av = avma;
    2608       19108 :   if (s < 0) x = RgXQ_inv(x, T);
    2609       19108 :   x = gen_pow_i(x, n, (void*)T, &_sqr, &_mul);
    2610       19108 :   return gerepilecopy(av, x);
    2611             : }
    2612             : static GEN
    2613      159050 : _ZXQsqr(void *data, GEN x) { return ZXQ_sqr(x, (GEN)data); }
    2614             : static GEN
    2615       50463 : _ZXQmul(void *data, GEN x, GEN y) { return ZXQ_mul(x,y, (GEN)data); }
    2616             : /* x,T in Z[X], n in N, compute lift(x^n mod T)) */
    2617             : GEN
    2618      163658 : ZXQ_powu(GEN x, ulong n, GEN T)
    2619             : {
    2620      163658 :   pari_sp av = avma;
    2621             : 
    2622      163658 :   if (!n) return pol_1(varn(x));
    2623      163658 :   if (n == 1) return ZX_copy(x);
    2624      107403 :   x = gen_powu_i(x, n, (void*)T, &_ZXQsqr, &_ZXQmul);
    2625      107399 :   return gerepilecopy(av, x);
    2626             : }
    2627             : 
    2628             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2629             : GEN
    2630        5306 : RgXQ_powers(GEN x, long l, GEN T)
    2631             : {
    2632        5306 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2633        5306 :   return gen_powers(x, l, use_sqr, (void *)T,_sqr,_mul,_one);
    2634             : }
    2635             : 
    2636             : /* a in K = Q[X]/(T), returns [a^0, ..., a^n] */
    2637             : GEN
    2638        3458 : QXQ_powers(GEN a, long n, GEN T)
    2639             : {
    2640        3458 :   GEN den, v = RgXQ_powers(Q_remove_denom(a, &den), n, T);
    2641             :   /* den*a integral; v[i+1] = (den*a)^i in K */
    2642        3458 :   if (den)
    2643             :   { /* restore denominators */
    2644        2135 :     GEN d = den;
    2645             :     long i;
    2646        2135 :     gel(v,2) = a;
    2647        7000 :     for (i=3; i<=n+1; i++) {
    2648        4865 :       d = mulii(d,den);
    2649        4865 :       gel(v,i) = RgX_Rg_div(gel(v,i), d);
    2650             :     }
    2651             :   }
    2652        3458 :   return v;
    2653             : }
    2654             : 
    2655             : static GEN
    2656        2016 : do_QXQ_eval(GEN v, long imin, GEN a, GEN T)
    2657             : {
    2658        2016 :   long l, i, m = 0;
    2659             :   GEN dz, z;
    2660        2016 :   GEN V = cgetg_copy(v, &l);
    2661        6496 :   for (i = imin; i < l; i++)
    2662             :   {
    2663        4480 :     GEN c = gel(v, i);
    2664        4480 :     if (typ(c) == t_POL) m = maxss(m, degpol(c));
    2665             :   }
    2666        2016 :   z = Q_remove_denom(QXQ_powers(a, m, T), &dz);
    2667        2079 :   for (i = 1; i < imin; i++) V[i] = v[i];
    2668        6496 :   for (i = imin; i < l; i++)
    2669             :   {
    2670        4480 :     GEN c = gel(v,i);
    2671        4480 :     if (typ(c) == t_POL) c = QX_ZXQV_eval(c, z, dz);
    2672        4480 :     gel(V,i) = c;
    2673             :   }
    2674        2016 :   return V;
    2675             : }
    2676             : /* [ s(a mod T) | s <- lift(v) ], a,T are QX, v a QXV */
    2677             : GEN
    2678        1953 : QXV_QXQ_eval(GEN v, GEN a, GEN T)
    2679        1953 : { return do_QXQ_eval(v, 1, a, T); }
    2680             : GEN
    2681          63 : QXX_QXQ_eval(GEN v, GEN a, GEN T)
    2682          63 : { return normalizepol(do_QXQ_eval(v, 2, a, T)); }
    2683             : 
    2684             : GEN
    2685         588 : RgXQ_matrix_pow(GEN y, long n, long m, GEN P)
    2686             : {
    2687         588 :   return RgXV_to_RgM(RgXQ_powers(y,m-1,P),n);
    2688             : }
    2689             : 
    2690             : GEN
    2691        5028 : RgXQ_norm(GEN x, GEN T)
    2692             : {
    2693             :   pari_sp av;
    2694        5028 :   long dx = degpol(x);
    2695             :   GEN L, y;
    2696             : 
    2697        5028 :   av = avma; y = resultant(T, x);
    2698        5028 :   L = leading_coeff(T);
    2699        5028 :   if (gequal1(L) || !signe(x)) return y;
    2700           0 :   return gerepileupto(av, gdiv(y, gpowgs(L, dx)));
    2701             : }
    2702             : 
    2703             : GEN
    2704     1622182 : RgX_blocks(GEN P, long n, long m)
    2705             : {
    2706     1622182 :   GEN z = cgetg(m+1,t_VEC);
    2707     1622182 :   long i,j, k=2, l = lg(P);
    2708     5059718 :   for(i=1; i<=m; i++)
    2709             :   {
    2710     3437536 :     GEN zi = cgetg(n+2,t_POL);
    2711     3437536 :     zi[1] = P[1];
    2712     3437536 :     gel(z,i) = zi;
    2713    11163125 :     for(j=2; j<n+2; j++)
    2714     7725589 :       gel(zi, j) = k==l ? gen_0 : gel(P,k++);
    2715     3437536 :     zi = RgX_renormalize_lg(zi, n+2);
    2716             :   }
    2717     1622182 :   return z;
    2718             : }
    2719             : 
    2720             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
    2721             : void
    2722     5895744 : RgX_even_odd(GEN p, GEN *pe, GEN *po)
    2723             : {
    2724     5895744 :   long n = degpol(p), v = varn(p), n0, n1, i;
    2725             :   GEN p0, p1;
    2726             : 
    2727     5895746 :   if (n <= 0) { *pe = RgX_copy(p); *po = zeropol(v); return; }
    2728             : 
    2729     5895459 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
    2730     5895459 :   p0 = cgetg(n0+2, t_POL); p0[1] = evalvarn(v)|evalsigne(1);
    2731     5895461 :   p1 = cgetg(n1+2, t_POL); p1[1] = evalvarn(v)|evalsigne(1);
    2732    24637241 :   for (i=0; i<n1; i++)
    2733             :   {
    2734    18741781 :     p0[2+i] = p[2+(i<<1)];
    2735    18741781 :     p1[2+i] = p[3+(i<<1)];
    2736             :   }
    2737     5895460 :   if (n1 != n0)
    2738     2881176 :     p0[2+i] = p[2+(i<<1)];
    2739     5895460 :   *pe = normalizepol(p0);
    2740     5895460 :   *po = normalizepol(p1);
    2741             : }
    2742             : 
    2743             : /* write p(X) = a_0(X^k) + Xa_1(X^k) + ... + X^(k-1)a_{k-1}(X^k), shallow function */
    2744             : GEN
    2745       40922 : RgX_splitting(GEN p, long k)
    2746             : {
    2747       40922 :   long n = degpol(p), v = varn(p), m, i, j, l;
    2748             :   GEN r;
    2749             : 
    2750       40922 :   m = n/k;
    2751       40922 :   r = cgetg(k+1,t_VEC);
    2752      225666 :   for(i=1; i<=k; i++)
    2753             :   {
    2754      184744 :     gel(r,i) = cgetg(m+3, t_POL);
    2755      184744 :     mael(r,i,1) = evalvarn(v)|evalsigne(1);
    2756             :   }
    2757      554981 :   for (j=1, i=0, l=2; i<=n; i++)
    2758             :   {
    2759      514059 :     gmael(r,j,l) = gel(p,2+i);
    2760      514059 :     if (j==k) { j=1; l++; } else j++;
    2761             :   }
    2762      225666 :   for(i=1; i<=k; i++)
    2763      184744 :     gel(r,i) = normalizepol_lg(gel(r,i),i<j?l+1:l);
    2764       40922 :   return r;
    2765             : }
    2766             : 
    2767             : /*******************************************************************/
    2768             : /*                                                                 */
    2769             : /*                        Kronecker form                           */
    2770             : /*                                                                 */
    2771             : /*******************************************************************/
    2772             : 
    2773             : /* z in R[Y] representing an elt in R[X,Y] mod T(Y) in Kronecker form,
    2774             :  * i.e subst(lift(z), x, y^(2deg(z)-1)). Recover the "real" z, with
    2775             :  * normalized coefficients */
    2776             : GEN
    2777      190997 : Kronecker_to_mod(GEN z, GEN T)
    2778             : {
    2779      190997 :   long i,j,lx,l = lg(z), N = (degpol(T)<<1) + 1;
    2780      190997 :   GEN x, t = cgetg(N,t_POL);
    2781      190997 :   t[1] = T[1];
    2782      190997 :   lx = (l-2) / (N-2); x = cgetg(lx+3,t_POL);
    2783      190997 :   x[1] = z[1];
    2784      190997 :   T = RgX_copy(T);
    2785     1880803 :   for (i=2; i<lx+2; i++, z+= N-2)
    2786             :   {
    2787     9205066 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
    2788     1689806 :     gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2789             :   }
    2790      190997 :   N = (l-2) % (N-2) + 2;
    2791      719034 :   for (j=2; j<N; j++) t[j] = z[j];
    2792      190997 :   gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2793      190997 :   return normalizepol_lg(x, i+1);
    2794             : }
    2795             : 
    2796             : /*******************************************************************/
    2797             : /*                                                                 */
    2798             : /*                        Domain detection                         */
    2799             : /*                                                                 */
    2800             : /*******************************************************************/
    2801             : 
    2802             : static GEN
    2803      151005 : zero_FpX_mod(GEN p, long v)
    2804             : {
    2805      151005 :   GEN r = cgetg(3,t_POL);
    2806      151005 :   r[1] = evalvarn(v);
    2807      151005 :   gel(r,2) = mkintmod(gen_0, icopy(p));
    2808      151005 :   return r;
    2809             : }
    2810             : 
    2811             : static GEN
    2812      499694 : RgX_mul_FpX(GEN x, GEN y, GEN p)
    2813             : {
    2814      499694 :   pari_sp av = avma;
    2815             :   GEN r;
    2816      499694 :   if (lgefint(p) == 3)
    2817             :   {
    2818      469562 :     ulong pp = uel(p, 2);
    2819      469562 :     r = Flx_to_ZX_inplace(Flx_mul(RgX_to_Flx(x, pp),
    2820             :                                   RgX_to_Flx(y, pp), pp));
    2821             :   }
    2822             :   else
    2823       30132 :     r = FpX_mul(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2824      499694 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2825      386377 :   return gerepileupto(av, FpX_to_mod(r, p));
    2826             : }
    2827             : 
    2828             : static GEN
    2829         532 : zero_FpXQX_mod(GEN pol, GEN p, long v)
    2830             : {
    2831         532 :   GEN r = cgetg(3,t_POL);
    2832         532 :   r[1] = evalvarn(v);
    2833         532 :   gel(r,2) = mkpolmod(mkintmod(gen_0, icopy(p)), gcopy(pol));
    2834         532 :   return r;
    2835             : }
    2836             : 
    2837             : static GEN
    2838        2037 : RgX_mul_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2839             : {
    2840        2037 :   pari_sp av = avma;
    2841             :   long dT;
    2842             :   GEN kx, ky, r;
    2843        2037 :   GEN T = RgX_to_FpX(pol, p);
    2844        2037 :   if (signe(T)==0) pari_err_OP("*", x, y);
    2845        2030 :   dT = degpol(T);
    2846        2030 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2847        2030 :   ky = ZXX_to_Kronecker(RgX_to_FpXQX(y, T, p), dT);
    2848        2030 :   r = FpX_mul(kx, ky, p);
    2849        2030 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2850        1505 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2851             : }
    2852             : 
    2853             : static GEN
    2854      514650 : RgX_liftred(GEN x, GEN T)
    2855      514650 : { return RgXQX_red(liftpol_shallow(x), T); }
    2856             : 
    2857             : static GEN
    2858      176989 : RgX_mul_QXQX(GEN x, GEN y, GEN T)
    2859             : {
    2860      176989 :   pari_sp av = avma;
    2861      176989 :   long dT = degpol(T);
    2862      176989 :   GEN r = QX_mul(ZXX_to_Kronecker(RgX_liftred(x, T), dT),
    2863             :                  ZXX_to_Kronecker(RgX_liftred(y, T), dT));
    2864      176989 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2865             : }
    2866             : 
    2867             : static GEN
    2868        1330 : RgX_sqr_FpX(GEN x, GEN p)
    2869             : {
    2870        1330 :   pari_sp av = avma;
    2871             :   GEN r;
    2872        1330 :   if (lgefint(p) == 3)
    2873             :   {
    2874        1126 :     ulong pp = uel(p, 2);
    2875        1126 :     r = Flx_to_ZX_inplace(Flx_sqr(RgX_to_Flx(x, pp), pp));
    2876             :   }
    2877             :   else
    2878         204 :     r = FpX_sqr(RgX_to_FpX(x, p), p);
    2879        1330 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2880        1232 :   return gerepileupto(av, FpX_to_mod(r, p));
    2881             : }
    2882             : 
    2883             : static GEN
    2884         196 : RgX_sqr_FpXQX(GEN x, GEN pol, GEN p)
    2885             : {
    2886         196 :   pari_sp av = avma;
    2887             :   long dT;
    2888         196 :   GEN kx, r, T = RgX_to_FpX(pol, p);
    2889         196 :   if (signe(T)==0) pari_err_OP("*",x,x);
    2890         189 :   dT = degpol(T);
    2891         189 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2892         189 :   r = FpX_sqr(kx, p);
    2893         189 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2894         189 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2895             : }
    2896             : 
    2897             : static GEN
    2898       12314 : RgX_sqr_QXQX(GEN x, GEN T)
    2899             : {
    2900       12314 :   pari_sp av = avma;
    2901       12314 :   long dT = degpol(T);
    2902       12314 :   GEN r = QX_sqr(ZXX_to_Kronecker(RgX_liftred(x, T), dT));
    2903       12314 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2904             : }
    2905             : 
    2906             : static GEN
    2907       67466 : RgX_rem_FpX(GEN x, GEN y, GEN p)
    2908             : {
    2909       67466 :   pari_sp av = avma;
    2910             :   GEN r;
    2911       67466 :   if (lgefint(p) == 3)
    2912             :   {
    2913       50676 :     ulong pp = uel(p, 2);
    2914       50676 :     r = Flx_to_ZX_inplace(Flx_rem(RgX_to_Flx(x, pp),
    2915             :                                   RgX_to_Flx(y, pp), pp));
    2916             :   }
    2917             :   else
    2918       16790 :     r = FpX_rem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2919       67466 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2920       29876 :   return gerepileupto(av, FpX_to_mod(r, p));
    2921             : }
    2922             : 
    2923             : static GEN
    2924       74179 : RgX_rem_QXQX(GEN x, GEN y, GEN T)
    2925             : {
    2926       74179 :   pari_sp av = avma;
    2927             :   GEN r;
    2928       74179 :   r = RgXQX_rem(RgX_liftred(x, T), RgX_liftred(y, T), T);
    2929       74179 :   return gerepilecopy(av, QXQX_to_mod_shallow(r, T));
    2930             : }
    2931             : static GEN
    2932          70 : RgX_rem_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2933             : {
    2934          70 :   pari_sp av = avma;
    2935             :   GEN r;
    2936          70 :   GEN T = RgX_to_FpX(pol, p);
    2937          70 :   if (signe(T) == 0) pari_err_OP("%", x, y);
    2938          63 :   if (lgefint(p) == 3)
    2939             :   {
    2940          55 :     ulong pp = uel(p, 2);
    2941          55 :     GEN Tp = ZX_to_Flx(T, pp);
    2942          55 :     r = FlxX_to_ZXX(FlxqX_rem(RgX_to_FlxqX(x, Tp, pp),
    2943             :                               RgX_to_FlxqX(y, Tp, pp), Tp, pp));
    2944             :   }
    2945             :   else
    2946           8 :     r = FpXQX_rem(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), T, p);
    2947          63 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2948          56 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2949             : }
    2950             : 
    2951             : #define code(t1,t2) ((t1 << 6) | t2)
    2952             : static GEN
    2953    53534429 : RgX_mul_fast(GEN x, GEN y)
    2954             : {
    2955             :   GEN p, pol;
    2956             :   long pa;
    2957    53534429 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2958    53534431 :   switch(t)
    2959             :   {
    2960    34543391 :     case t_INT:    return ZX_mul(x,y);
    2961     2654214 :     case t_FRAC:   return QX_mul(x,y);
    2962      104693 :     case t_FFELT:  return FFX_mul(x, y, pol);
    2963      499694 :     case t_INTMOD: return RgX_mul_FpX(x, y, p);
    2964      176989 :     case code(t_POLMOD, t_INT):
    2965             :     case code(t_POLMOD, t_FRAC):
    2966      176989 :                    return RgX_mul_QXQX(x, y, pol);
    2967        2037 :     case code(t_POLMOD, t_INTMOD):
    2968        2037 :                    return RgX_mul_FpXQX(x, y, pol, p);
    2969    15553413 :     default:       return NULL;
    2970             :   }
    2971             : }
    2972             : static GEN
    2973      139262 : RgX_sqr_fast(GEN x)
    2974             : {
    2975             :   GEN p, pol;
    2976             :   long pa;
    2977      139262 :   long t = RgX_type(x,&p,&pol,&pa);
    2978      139262 :   switch(t)
    2979             :   {
    2980      103806 :     case t_INT:    return ZX_sqr(x);
    2981       11410 :     case t_FRAC:   return QX_sqr(x);
    2982        2296 :     case t_FFELT:  return FFX_sqr(x, pol);
    2983        1330 :     case t_INTMOD: return RgX_sqr_FpX(x, p);
    2984       12314 :     case code(t_POLMOD, t_INT):
    2985             :     case code(t_POLMOD, t_FRAC):
    2986       12314 :                    return RgX_sqr_QXQX(x, pol);
    2987         196 :     case code(t_POLMOD, t_INTMOD):
    2988         196 :                    return RgX_sqr_FpXQX(x, pol, p);
    2989        7910 :     default:       return NULL;
    2990             :   }
    2991             : }
    2992             : 
    2993             : static GEN
    2994     7530569 : RgX_rem_fast(GEN x, GEN y)
    2995             : {
    2996             :   GEN p, pol;
    2997             :   long pa;
    2998     7530569 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2999     7530569 :   switch(t)
    3000             :   {
    3001     4478834 :     case t_INT:    return ZX_is_monic(y) ? ZX_rem(x,y): NULL;
    3002     1603156 :     case t_FRAC:   return RgX_is_ZX(y) && ZX_is_monic(y) ? QX_ZX_rem(x,y): NULL;
    3003          84 :     case t_FFELT:  return FFX_rem(x, y, pol);
    3004       67466 :     case t_INTMOD: return RgX_rem_FpX(x, y, p);
    3005       74179 :     case code(t_POLMOD, t_INT):
    3006             :     case code(t_POLMOD, t_FRAC):
    3007       74179 :                    return RgX_rem_QXQX(x, y, pol);
    3008          70 :     case code(t_POLMOD, t_INTMOD):
    3009          70 :                    return RgX_rem_FpXQX(x, y, pol, p);
    3010     1306780 :     default:       return NULL;
    3011             :   }
    3012             : }
    3013             : 
    3014             : #undef code
    3015             : 
    3016             : GEN
    3017    45641243 : RgX_mul(GEN x, GEN y)
    3018             : {
    3019    45641243 :   GEN z = RgX_mul_fast(x,y);
    3020    45641236 :   if (!z) z = RgX_mul_i(x,y);
    3021    45641236 :   return z;
    3022             : }
    3023             : 
    3024             : GEN
    3025      127404 : RgX_sqr(GEN x)
    3026             : {
    3027      127404 :   GEN z = RgX_sqr_fast(x);
    3028      127397 :   if (!z) z = RgX_sqr_i(x);
    3029      127397 :   return z;
    3030             : }
    3031             : 
    3032             : GEN
    3033     7530569 : RgX_rem(GEN x, GEN y)
    3034             : {
    3035     7530569 :   GEN z = RgX_rem_fast(x, y);
    3036     7530562 :   if (!z) z = RgX_divrem_i(x, y, ONLY_REM);
    3037     7530562 :   return z;
    3038             : }
    3039             : 
    3040             : GEN
    3041     6403844 : RgXn_mul(GEN f, GEN g, long n)
    3042             : {
    3043     6403844 :   pari_sp av = avma;
    3044     6403844 :   GEN h = RgX_mul_fast(f,g);
    3045     6403844 :   if (!h) return RgXn_mul2(f,g,n);
    3046     1251571 :   if (degpol(h) < n) return h;
    3047      480472 :   return gerepilecopy(av, RgXn_red_shallow(h, n));
    3048             : }
    3049             : 
    3050             : GEN
    3051       11858 : RgXn_sqr(GEN f, long n)
    3052             : {
    3053       11858 :   pari_sp av = avma;
    3054       11858 :   GEN g = RgX_sqr_fast(f);
    3055       11858 :   if (!g) return RgXn_sqr2(f,n);
    3056       11613 :   if (degpol(g) < n) return g;
    3057       10521 :   return gerepilecopy(av, RgXn_red_shallow(g, n));
    3058             : }
    3059             : 
    3060             : /* (f*g) \/ x^n */
    3061             : GEN
    3062     1489342 : RgX_mulhigh_i(GEN f, GEN g, long n)
    3063             : {
    3064     1489342 :   GEN h = RgX_mul_fast(f,g);
    3065     1489342 :   return h? RgX_shift_shallow(h, -n): RgX_mulhigh_i2(f,g,n);
    3066             : }
    3067             : 
    3068             : /* (f*g) \/ x^n */
    3069             : GEN
    3070           0 : RgX_sqrhigh_i(GEN f, long n)
    3071             : {
    3072           0 :   GEN h = RgX_sqr_fast(f);
    3073           0 :   return h? RgX_shift_shallow(h, -n): RgX_sqrhigh_i2(f,n);
    3074             : }

Generated by: LCOV version 1.13