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 25579-8c4672f557) Lines: 1531 1715 89.3 %
Date: 2020-07-09 06:03:45 Functions: 180 199 90.5 %
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     4232412 : brent_kung_optpow(long d, long n, long m)
      29             : {
      30             :   long p, r;
      31     4232412 :   long pold=1, rold=n*(d-1);
      32    26428160 :   for(p=2; p<=d; p++)
      33             :   {
      34    22195748 :     r = m*(p-1) + n*((d-1)/p);
      35    22195748 :     if (r<rold) { pold=p; rold=r; }
      36             :   }
      37     4232412 :   return pold;
      38             : }
      39             : 
      40             : static GEN
      41     7131289 : 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     7131289 :   pari_sp av = avma;
      45             :   long i;
      46     7131289 :   GEN z = cmul(E,P,a,ff->one(E));
      47     7130729 :   if (!z) z = gen_0;
      48    51182864 :   for (i=1; i<=n; i++)
      49             :   {
      50    44052435 :     GEN t = cmul(E,P,a+i,gel(V,i+1));
      51    44061458 :     if (t) {
      52    32839090 :       z = ff->add(E, z, t);
      53    32816375 :       if (gc_needed(av,2)) z = gerepileupto(av, z);
      54             :     }
      55             :   }
      56     7130429 :   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     5922979 : 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     5922979 :   pari_sp av = avma;
      69     5922979 :   long l = lg(V)-1;
      70             :   GEN z, u;
      71             : 
      72     5922979 :   if (d < 0) return ff->zero(E);
      73     5604975 :   if (d < l) return gerepileupto(av, gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul));
      74      656602 :   if (l<2) pari_err_DOMAIN("gen_RgX_bkeval_powers", "#powers", "<",gen_2,V);
      75      656602 :   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      656602 :   d -= l;
      81      656602 :   z = gen_RgXQ_eval_powers(P,V,d+1,l-1,E,ff,cmul);
      82     1526816 :   while (d >= l-1)
      83             :   {
      84      870243 :     d -= l-1;
      85      870243 :     u = gen_RgXQ_eval_powers(P,V,d+1,l-2,E,ff,cmul);
      86      870096 :     z = ff->add(E,u, ff->mul(E,z,gel(V,l)));
      87      870161 :     if (gc_needed(av,2))
      88          91 :       z = gerepileupto(av, z);
      89             :   }
      90      656573 :   u = gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul);
      91      656593 :   z = ff->add(E,u, ff->mul(E,z,gel(V,d+2)));
      92      656579 :   return gerepileupto(av, ff->red(E,z));
      93             : }
      94             : 
      95             : GEN
      96      425829 : 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      425829 :   pari_sp av = avma;
     100             :   GEN z, V;
     101             :   long rtd;
     102      425829 :   if (d < 0) return ff->zero(E);
     103      425255 :   rtd = (long) sqrt((double)d);
     104      425255 :   V = gen_powers(x,rtd,use_sqr,E,ff->sqr,ff->mul,ff->one);
     105      425259 :   z = gen_bkeval_powers(Q, d, V, E, ff, cmul);
     106      425266 :   return gerepileupto(av, z);
     107             : }
     108             : 
     109             : static GEN
     110     1092368 : _gen_nored(void *E, GEN x) { (void)E; return x; }
     111             : static GEN
     112    13978187 : _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      941367 : _gen_mul(void *E, GEN x, GEN y) { (void)E; return gmul(x, y); }
     117             : static GEN
     118      308294 : _gen_sqr(void *E, GEN x) { (void)E; return gsqr(x); }
     119             : static GEN
     120     1114966 : _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 : QXQX_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 :     GEN c = gel(P,i+2), b = gel(B,d+1-i);
     190        3661 :     s = RgX_add(QXQX_mul(s, A, T), typ(c)==t_POL ? QXQX_QXQ_mul(b, c, T): gmul(b, c));
     191        3661 :     if (gc_needed(av,1))
     192             :     {
     193           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"QXQX_homogenous_eval(%ld)",i);
     194           0 :       s = gerepileupto(av, s);
     195             :     }
     196             :   }
     197        1568 :   s = gerepileupto(av, s);
     198        1568 :   return mkvec2(s, gel(B,d+1));
     199             : }
     200             : 
     201             : const struct bb_algebra *
     202      142984 : get_Rg_algebra(void)
     203             : {
     204      142984 :   return &Rg_algebra;
     205             : }
     206             : 
     207             : static struct bb_ring Rg_ring = {  _gen_add, _gen_mul, _gen_sqr };
     208             : 
     209             : static GEN
     210       17780 : _RgX_divrem(void *E, GEN x, GEN y, GEN *r)
     211             : {
     212             :   (void) E;
     213       17780 :   return RgX_divrem(x, y, r);
     214             : }
     215             : 
     216             : GEN
     217        4382 : RgX_digits(GEN x, GEN T)
     218             : {
     219        4382 :   pari_sp av = avma;
     220        4382 :   long d = degpol(T), n = (lgpol(x)+d-1)/d;
     221        4382 :   GEN z = gen_digits(x,T,n,NULL, &Rg_ring, _RgX_divrem);
     222        4382 :   return gerepileupto(av, z);
     223             : }
     224             : 
     225             : /*******************************************************************/
     226             : /*                                                                 */
     227             : /*                         RgX                                     */
     228             : /*                                                                 */
     229             : /*******************************************************************/
     230             : 
     231             : long
     232    22180873 : RgX_equal(GEN x, GEN y)
     233             : {
     234    22180873 :   long i = lg(x);
     235             : 
     236    22180873 :   if (i != lg(y)) return 0;
     237    94240124 :   for (i--; i > 1; i--)
     238    72188946 :     if (!gequal(gel(x,i),gel(y,i))) return 0;
     239    22051178 :   return 1;
     240             : }
     241             : 
     242             : /* Returns 1 in the base ring over which x is defined */
     243             : /* HACK: this also works for t_SER */
     244             : GEN
     245     3341045 : Rg_get_1(GEN x)
     246             : {
     247             :   GEN p, T;
     248     3341045 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     249     3341045 :   if (RgX_type_is_composite(tx))
     250       32725 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     251     3341045 :   switch(tx)
     252             :   {
     253         357 :     case t_INTMOD: retmkintmod(is_pm1(p)? gen_0: gen_1, icopy(p));
     254           7 :     case t_PADIC: return cvtop(gen_1, p, lx);
     255         728 :     case t_FFELT: return FF_1(T);
     256     3339953 :     default: return gen_1;
     257             :   }
     258             : }
     259             : /* Returns 0 in the base ring over which x is defined */
     260             : /* HACK: this also works for t_SER */
     261             : GEN
     262     6972441 : Rg_get_0(GEN x)
     263             : {
     264             :   GEN p, T;
     265     6972441 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     266     6972441 :   if (RgX_type_is_composite(tx))
     267       51310 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     268     6972441 :   switch(tx)
     269             :   {
     270         441 :     case t_INTMOD: retmkintmod(gen_0, icopy(p));
     271           0 :     case t_PADIC: return zeropadic(p, lx);
     272         266 :     case t_FFELT: return FF_zero(T);
     273     6971734 :     default: return gen_0;
     274             :   }
     275             : }
     276             : 
     277             : GEN
     278        4298 : QX_ZXQV_eval(GEN P, GEN V, GEN dV)
     279             : {
     280        4298 :   long i, n = degpol(P);
     281             :   GEN z, dz, dP;
     282        4298 :   if (n < 0) return gen_0;
     283        4298 :   P = Q_remove_denom(P, &dP);
     284        4298 :   z = gel(P,2); if (n == 0) return icopy(z);
     285        2345 :   if (dV) z = mulii(dV, z); /* V[1] = dV */
     286        2345 :   z = ZX_Z_add_shallow(ZX_Z_mul(gel(V,2),gel(P,3)), z);
     287        4949 :   for (i=2; i<=n; i++) z = ZX_add(ZX_Z_mul(gel(V,i+1),gel(P,2+i)), z);
     288        2345 :   dz = mul_denom(dP, dV);
     289        2345 :   return dz? RgX_Rg_div(z, dz): z;
     290             : }
     291             : 
     292             : /* Return P(h * x), not memory clean */
     293             : GEN
     294        6930 : RgX_unscale(GEN P, GEN h)
     295             : {
     296        6930 :   long i, l = lg(P);
     297        6930 :   GEN hi = gen_1, Q = cgetg(l, t_POL);
     298        6930 :   Q[1] = P[1];
     299        6930 :   if (l == 2) return Q;
     300        6902 :   gel(Q,2) = gcopy(gel(P,2));
     301       37184 :   for (i=3; i<l; i++)
     302             :   {
     303       30282 :     hi = gmul(hi,h);
     304       30282 :     gel(Q,i) = gmul(gel(P,i), hi);
     305             :   }
     306        6902 :   return Q;
     307             : }
     308             : /* P a ZX, Return P(h * x), not memory clean; optimize for h = -1 */
     309             : GEN
     310      373284 : ZX_z_unscale(GEN P, long h)
     311             : {
     312      373284 :   long i, l = lg(P);
     313      373284 :   GEN Q = cgetg(l, t_POL);
     314      373284 :   Q[1] = P[1];
     315      373284 :   if (l == 2) return Q;
     316      373284 :   gel(Q,2) = gel(P,2);
     317      373284 :   if (l == 3) return Q;
     318      373284 :   if (h == -1)
     319      324856 :     for (i = 3; i < l; i++)
     320             :     {
     321      312862 :       gel(Q,i) = negi(gel(P,i));
     322      312862 :       if (++i == l) break;
     323      309041 :       gel(Q,i) = gel(P,i);
     324             :     }
     325             :   else
     326             :   {
     327             :     GEN hi;
     328      357469 :     gel(Q,3) = mulis(gel(P,3), h);
     329      357469 :     hi = sqrs(h);
     330     1724016 :     for (i = 4; i < l; i++)
     331             :     {
     332     1366547 :       gel(Q,i) = mulii(gel(P,i), hi);
     333     1366547 :       if (i != l-1) hi = mulis(hi,h);
     334             :     }
     335             :   }
     336      373284 :   return Q;
     337             : }
     338             : /* P a ZX, h a t_INT. Return P(h * x), not memory clean; optimize for h = -1 */
     339             : GEN
     340        7273 : ZX_unscale(GEN P, GEN h)
     341             : {
     342             :   long i, l;
     343             :   GEN Q, hi;
     344        7273 :   i = itos_or_0(h); if (i) return ZX_z_unscale(P, i);
     345           7 :   l = lg(P); Q = cgetg(l, t_POL);
     346           7 :   Q[1] = P[1];
     347           7 :   if (l == 2) return Q;
     348           7 :   gel(Q,2) = gel(P,2);
     349           7 :   if (l == 3) return Q;
     350           7 :   hi = h;
     351           7 :   gel(Q,3) = mulii(gel(P,3), hi);
     352          70 :   for (i = 4; i < l; i++)
     353             :   {
     354          63 :     hi = mulii(hi,h);
     355          63 :     gel(Q,i) = mulii(gel(P,i), hi);
     356             :   }
     357           7 :   return Q;
     358             : }
     359             : /* P a ZX. Return P(x << n), not memory clean */
     360             : GEN
     361       64483 : ZX_unscale2n(GEN P, long n)
     362             : {
     363       64483 :   long i, ni = n, l = lg(P);
     364       64483 :   GEN Q = cgetg(l, t_POL);
     365       64483 :   Q[1] = P[1];
     366       64483 :   if (l == 2) return Q;
     367       64483 :   gel(Q,2) = gel(P,2);
     368       64483 :   if (l == 3) return Q;
     369       64483 :   gel(Q,3) = shifti(gel(P,3), ni);
     370      258082 :   for (i=4; i<l; i++)
     371             :   {
     372      193599 :     ni += n;
     373      193599 :     gel(Q,i) = shifti(gel(P,i), ni);
     374             :   }
     375       64483 :   return Q;
     376             : }
     377             : /* P(h*X) / h, assuming h | P(0), i.e. the result is a ZX */
     378             : GEN
     379        1190 : ZX_unscale_div(GEN P, GEN h)
     380             : {
     381        1190 :   long i, l = lg(P);
     382        1190 :   GEN hi, Q = cgetg(l, t_POL);
     383        1190 :   Q[1] = P[1];
     384        1190 :   if (l == 2) return Q;
     385        1190 :   gel(Q,2) = diviiexact(gel(P,2), h);
     386        1190 :   if (l == 3) return Q;
     387        1190 :   gel(Q,3) = gel(P,3);
     388        1190 :   if (l == 4) return Q;
     389        1190 :   hi = h;
     390        1190 :   gel(Q,4) = mulii(gel(P,4), hi);
     391        5257 :   for (i=5; i<l; i++)
     392             :   {
     393        4067 :     hi = mulii(hi,h);
     394        4067 :     gel(Q,i) = mulii(gel(P,i), hi);
     395             :   }
     396        1190 :   return Q;
     397             : }
     398             : 
     399             : GEN
     400         658 : RgXV_unscale(GEN v, GEN h)
     401             : {
     402             :   long i, l;
     403             :   GEN w;
     404         658 :   if (!h || isint1(h)) return v;
     405         350 :   w = cgetg_copy(v, &l);
     406        1435 :   for (i=1; i<l; i++) gel(w,i) = RgX_unscale(gel(v,i), h);
     407         350 :   return w;
     408             : }
     409             : 
     410             : /* Return h^degpol(P) P(x / h), not memory clean */
     411             : GEN
     412      513359 : RgX_rescale(GEN P, GEN h)
     413             : {
     414      513359 :   long i, l = lg(P);
     415      513359 :   GEN Q = cgetg(l,t_POL), hi = h;
     416      513359 :   gel(Q,l-1) = gel(P,l-1);
     417     1871714 :   for (i=l-2; i>=2; i--)
     418             :   {
     419     1869425 :     gel(Q,i) = gmul(gel(P,i), hi);
     420     1869425 :     if (i == 2) break;
     421     1358355 :     hi = gmul(hi,h);
     422             :   }
     423      513359 :   Q[1] = P[1]; return Q;
     424             : }
     425             : 
     426             : /* A(X^d) --> A(X) */
     427             : GEN
     428      160256 : RgX_deflate(GEN x0, long d)
     429             : {
     430             :   GEN z, y, x;
     431      160256 :   long i,id, dy, dx = degpol(x0);
     432      160256 :   if (d == 1 || dx <= 0) return leafcopy(x0);
     433       84501 :   dy = dx/d;
     434       84501 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     435       84501 :   z = y + 2;
     436       84501 :   x = x0+ 2;
     437      375900 :   for (i=id=0; i<=dy; i++,id+=d) gel(z,i) = gel(x,id);
     438       84501 :   return y;
     439             : }
     440             : 
     441             : /* F a t_RFRAC */
     442             : long
     443          70 : rfrac_deflate_order(GEN F)
     444             : {
     445          70 :   GEN N = gel(F,1), D = gel(F,2);
     446          70 :   long m = (degpol(D) <= 0)? 0: RgX_deflate_order(D);
     447          70 :   if (m == 1) return 1;
     448          28 :   if (typ(N) == t_POL && varn(N) == varn(D))
     449          14 :     m = cgcd(m, RgX_deflate_order(N));
     450          28 :   return m;
     451             : }
     452             : /* F a t_RFRAC */
     453             : GEN
     454          70 : rfrac_deflate_max(GEN F, long *m)
     455             : {
     456          70 :   *m = rfrac_deflate_order(F);
     457          70 :   return rfrac_deflate(F, *m);
     458             : }
     459             : /* F a t_RFRAC */
     460             : GEN
     461          70 : rfrac_deflate(GEN F, long m)
     462             : {
     463          70 :   GEN N = gel(F,1), D = gel(F,2);
     464          70 :   if (m == 1) return F;
     465          28 :   if (typ(N) == t_POL && varn(N) == varn(D)) N = RgX_deflate(N, m);
     466          28 :   D = RgX_deflate(D, m); return mkrfrac(N, D);
     467             : }
     468             : 
     469             : /* return x0(X^d) */
     470             : GEN
     471      601566 : RgX_inflate(GEN x0, long d)
     472             : {
     473      601566 :   long i, id, dy, dx = degpol(x0);
     474      601566 :   GEN x = x0 + 2, z, y;
     475      601566 :   if (dx <= 0) return leafcopy(x0);
     476      593901 :   dy = dx*d;
     477      593901 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     478      593901 :   z = y + 2;
     479    21664156 :   for (i=0; i<=dy; i++) gel(z,i) = gen_0;
     480     9730761 :   for (i=id=0; i<=dx; i++,id+=d) gel(z,id) = gel(x,i);
     481      593901 :   return y;
     482             : }
     483             : 
     484             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     485             : static GEN
     486     1719250 : RgX_translate_basecase(GEN P, GEN c)
     487             : {
     488     1719250 :   pari_sp av = avma;
     489             :   GEN Q, R;
     490             :   long i, k, n;
     491             : 
     492     1719250 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     493     1716832 :   Q = leafcopy(P);
     494     1716832 :   R = Q+2; n = degpol(P);
     495     1716832 :   if (isint1(c))
     496             :   {
     497       10358 :     for (i=1; i<=n; i++)
     498             :     {
     499       36713 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gel(R,k+1));
     500        7846 :       if (gc_needed(av,2))
     501             :       {
     502           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(1), i = %ld/%ld", i,n);
     503           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     504             :       }
     505             :     }
     506             :   }
     507     1714320 :   else if (isintm1(c))
     508             :   {
     509      173726 :     for (i=1; i<=n; i++)
     510             :     {
     511     1124265 :       for (k=n-i; k<n; k++) gel(R,k) = gsub(gel(R,k), gel(R,k+1));
     512      148586 :       if (gc_needed(av,2))
     513             :       {
     514           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(-1), i = %ld/%ld", i,n);
     515           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     516             :       }
     517             :     }
     518             :   }
     519             :   else
     520             :   {
     521     5655500 :     for (i=1; i<=n; i++)
     522             :     {
     523    13274044 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gmul(c, gel(R,k+1)));
     524     3966320 :       if (gc_needed(av,2))
     525             :       {
     526           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate, i = %ld/%ld", i,n);
     527           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     528             :       }
     529             :     }
     530             :   }
     531     1716832 :   return gerepilecopy(av, Q);
     532             : }
     533             : GEN
     534     1719572 : RgX_translate(GEN P, GEN c)
     535             : {
     536     1719572 :   pari_sp av = avma;
     537     1719572 :   long n = degpol(P);
     538     1719572 :   if (n < 40)
     539     1719250 :     return RgX_translate_basecase(P, c);
     540             :   else
     541             :   {
     542         322 :     long d = n >> 1;
     543         322 :     GEN Q = RgX_translate(RgX_shift_shallow(P, -d), c);
     544         322 :     GEN R = RgX_translate(RgXn_red_shallow(P, d), c);
     545         322 :     GEN S = gpowgs(deg1pol_shallow(gen_1, c, varn(P)), d);
     546         322 :     return gerepileupto(av, RgX_add(RgX_mul(Q, S), R));
     547             :   }
     548             : }
     549             : 
     550             : /* return lift( P(X + c) ) using Horner, c in R[y]/(T) */
     551             : GEN
     552       21280 : RgXQX_translate(GEN P, GEN c, GEN T)
     553             : {
     554       21280 :   pari_sp av = avma;
     555             :   GEN Q, R;
     556             :   long i, k, n;
     557             : 
     558       21280 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     559       20951 :   Q = leafcopy(P);
     560       20951 :   R = Q+2; n = degpol(P);
     561       72842 :   for (i=1; i<=n; i++)
     562             :   {
     563      218806 :     for (k=n-i; k<n; k++)
     564             :     {
     565      166915 :       pari_sp av2 = avma;
     566      166915 :       gel(R,k) = gerepileupto(av2,
     567      166915 :                    RgX_rem(gadd(gel(R,k), gmul(c, gel(R,k+1))), T));
     568             :     }
     569       51891 :     if (gc_needed(av,2))
     570             :     {
     571           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXQX_translate, i = %ld/%ld", i,n);
     572           0 :       Q = gerepilecopy(av, Q); R = Q+2;
     573             :     }
     574             :   }
     575       20951 :   return gerepilecopy(av, Q);
     576             : }
     577             : 
     578             : /********************************************************************/
     579             : /**                                                                **/
     580             : /**                          CONVERSIONS                           **/
     581             : /**                       (not memory clean)                       **/
     582             : /**                                                                **/
     583             : /********************************************************************/
     584             : /* to INT / FRAC / (POLMOD mod T), not memory clean because T not copied,
     585             :  * but everything else is */
     586             : static GEN
     587       41860 : QXQ_to_mod(GEN x, GEN T)
     588             : {
     589             :   long d;
     590       41860 :   switch(typ(x))
     591             :   {
     592       17717 :     case t_INT:  return icopy(x);
     593         805 :     case t_FRAC: return gcopy(x);
     594       23338 :     case t_POL:
     595       23338 :       d = degpol(x);
     596       23338 :       if (d < 0) return gen_0;
     597       22470 :       if (d == 0) return gcopy(gel(x,2));
     598       21798 :       return mkpolmod(RgX_copy(x), T);
     599           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     600             :              return NULL;/* LCOV_EXCL_LINE */
     601             :   }
     602             : }
     603             : /* pure shallow version */
     604             : GEN
     605     1150360 : QXQ_to_mod_shallow(GEN x, GEN T)
     606             : {
     607             :   long d;
     608     1150360 :   switch(typ(x))
     609             :   {
     610      438242 :     case t_INT:
     611      438242 :     case t_FRAC: return x;
     612      712118 :     case t_POL:
     613      712118 :       d = degpol(x);
     614      712118 :       if (d < 0) return gen_0;
     615      538019 :       if (d == 0) return gel(x,2);
     616      515877 :       return mkpolmod(x, T);
     617           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     618             :              return NULL;/* LCOV_EXCL_LINE */
     619             :   }
     620             : }
     621             : /* T a ZX, z lifted from (Q[Y]/(T(Y)))[X], apply QXQ_to_mod to all coeffs.
     622             :  * Not memory clean because T not copied, but everything else is */
     623             : static GEN
     624        8876 : QXQX_to_mod(GEN z, GEN T)
     625             : {
     626        8876 :   long i,l = lg(z);
     627        8876 :   GEN x = cgetg(l,t_POL);
     628       44142 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod(gel(z,i), T);
     629        8876 :   x[1] = z[1]; return normalizepol_lg(x,l);
     630             : }
     631             : /* pure shallow version */
     632             : GEN
     633      150703 : QXQX_to_mod_shallow(GEN z, GEN T)
     634             : {
     635      150703 :   long i,l = lg(z);
     636      150703 :   GEN x = cgetg(l,t_POL);
     637      750758 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod_shallow(gel(z,i), T);
     638      150703 :   x[1] = z[1]; return normalizepol_lg(x,l);
     639             : }
     640             : /* Apply QXQX_to_mod to all entries. Memory-clean ! */
     641             : GEN
     642        2695 : QXQXV_to_mod(GEN V, GEN T)
     643             : {
     644        2695 :   long i, l = lg(V);
     645        2695 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     646       11571 :   for (i=1;i<l; i++) gel(z,i) = QXQX_to_mod(gel(V,i), T);
     647        2695 :   return z;
     648             : }
     649             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     650             : GEN
     651        5131 : QXQV_to_mod(GEN V, GEN T)
     652             : {
     653        5131 :   long i, l = lg(V);
     654        5131 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     655       11725 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod(gel(V,i), T);
     656        5131 :   return z;
     657             : }
     658             : 
     659             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     660             : GEN
     661       24465 : QXQC_to_mod_shallow(GEN V, GEN T)
     662             : {
     663       24465 :   long i, l = lg(V);
     664       24465 :   GEN z = cgetg(l, t_COL);
     665      574770 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod_shallow(gel(V,i), T);
     666       24465 :   return z;
     667             : }
     668             : 
     669             : GEN
     670        7476 : QXQM_to_mod_shallow(GEN V, GEN T)
     671             : {
     672        7476 :   long i, l = lg(V);
     673        7476 :   GEN z = cgetg(l, t_MAT);
     674       31941 :   for (i=1; i<l; i++) gel(z,i) = QXQC_to_mod_shallow(gel(V,i), T);
     675        7476 :   return z;
     676             : }
     677             : 
     678             : GEN
     679     3827012 : RgX_renormalize_lg(GEN x, long lx)
     680             : {
     681             :   long i;
     682     6383630 :   for (i = lx-1; i>1; i--)
     683     5997981 :     if (! gequal0(gel(x,i))) break; /* _not_ isexactzero */
     684     3827012 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
     685     3827012 :   setlg(x, i+1); setsigne(x, i != 1); return x;
     686             : }
     687             : 
     688             : GEN
     689      538872 : RgV_to_RgX(GEN x, long v)
     690             : {
     691      538872 :   long i, k = lg(x);
     692             :   GEN p;
     693             : 
     694     1963055 :   while (--k && gequal0(gel(x,k)));
     695      538872 :   if (!k) return pol_0(v);
     696      532327 :   i = k+2; p = cgetg(i,t_POL);
     697      532327 :   p[1] = evalsigne(1) | evalvarn(v);
     698     5672874 :   x--; for (k=2; k<i; k++) gel(p,k) = gel(x,k);
     699      532327 :   return p;
     700             : }
     701             : GEN
     702      159428 : RgV_to_RgX_reverse(GEN x, long v)
     703             : {
     704      159428 :   long j, k, l = lg(x);
     705             :   GEN p;
     706             : 
     707      159428 :   for (k = 1; k < l; k++)
     708      159428 :     if (!gequal0(gel(x,k))) break;
     709      159428 :   if (k == l) return pol_0(v);
     710      159428 :   k -= 1;
     711      159428 :   l -= k;
     712      159428 :   x += k;
     713      159428 :   p = cgetg(l+1,t_POL);
     714      159428 :   p[1] = evalsigne(1) | evalvarn(v);
     715      910437 :   for (j=2, k=l; j<=l; j++) gel(p,j) = gel(x,--k);
     716      159428 :   return p;
     717             : }
     718             : 
     719             : /* return the (N-dimensional) vector of coeffs of p */
     720             : GEN
     721     8478373 : RgX_to_RgC(GEN x, long N)
     722             : {
     723             :   long i, l;
     724             :   GEN z;
     725     8478373 :   l = lg(x)-1; x++;
     726     8478373 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
     727     8478373 :   z = cgetg(N+1,t_COL);
     728    54478304 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
     729    16683342 :   for (   ; i<=N; i++) gel(z,i) = gen_0;
     730     8478419 :   return z;
     731             : }
     732             : GEN
     733      748208 : Rg_to_RgC(GEN x, long N)
     734             : {
     735      748208 :   return (typ(x) == t_POL)? RgX_to_RgC(x,N): scalarcol_shallow(x, N);
     736             : }
     737             : 
     738             : /* vector of polynomials (in v) whose coeffs are given by the columns of x */
     739             : GEN
     740       62624 : RgM_to_RgXV(GEN x, long v)
     741      350784 : { pari_APPLY_type(t_VEC, RgV_to_RgX(gel(x,i), v)) }
     742             : 
     743             : /* matrix whose entries are given by the coeffs of the polynomials in
     744             :  * vector v (considered as degree n-1 polynomials) */
     745             : GEN
     746      131815 : RgV_to_RgM(GEN x, long n)
     747      876371 : { pari_APPLY_type(t_MAT, Rg_to_RgC(gel(x,i), n)) }
     748             : 
     749             : GEN
     750        6098 : RgXV_to_RgM(GEN x, long n)
     751       38420 : { pari_APPLY_type(t_MAT, RgX_to_RgC(gel(x,i), n)) }
     752             : 
     753             : /* polynomial (in v) of polynomials (in w) whose coeffs are given by the columns of x */
     754             : GEN
     755       17663 : RgM_to_RgXX(GEN x, long v,long w)
     756             : {
     757       17663 :   long j, lx = lg(x);
     758       17663 :   GEN y = cgetg(lx+1, t_POL);
     759       17663 :   y[1] = evalsigne(1) | evalvarn(v);
     760       17663 :   y++;
     761      100630 :   for (j=1; j<lx; j++) gel(y,j) = RgV_to_RgX(gel(x,j), w);
     762       17663 :   return normalizepol_lg(--y, lx+1);
     763             : }
     764             : 
     765             : /* matrix whose entries are given by the coeffs of the polynomial v in
     766             :  * two variables (considered as degree n-1 polynomials) */
     767             : GEN
     768         308 : RgXX_to_RgM(GEN v, long n)
     769             : {
     770         308 :   long j, N = lg(v)-1;
     771         308 :   GEN y = cgetg(N, t_MAT);
     772        1001 :   for (j=1; j<N; j++) gel(y,j) = Rg_to_RgC(gel(v,j+1), n);
     773         308 :   return y;
     774             : }
     775             : 
     776             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     777             : GEN
     778       21810 : RgXY_swapspec(GEN x, long n, long w, long nx)
     779             : {
     780       21810 :   long j, ly = n+3;
     781       21810 :   GEN y = cgetg(ly, t_POL);
     782       21810 :   y[1] = evalsigne(1);
     783      267543 :   for (j=2; j<ly; j++)
     784             :   {
     785             :     long k;
     786      245733 :     GEN a = cgetg(nx+2,t_POL);
     787      245733 :     a[1] = evalsigne(1) | evalvarn(w);
     788     1238754 :     for (k=0; k<nx; k++)
     789             :     {
     790      993021 :       GEN xk = gel(x,k);
     791      993021 :       if (typ(xk)==t_POL)
     792      914086 :         gel(a,k+2) = j<lg(xk)? gel(xk,j): gen_0;
     793             :       else
     794       78935 :         gel(a,k+2) = j==2 ? xk: gen_0;
     795             :     }
     796      245733 :     gel(y,j) = normalizepol_lg(a, nx+2);
     797             :   }
     798       21810 :   return normalizepol_lg(y,ly);
     799             : }
     800             : 
     801             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     802             : GEN
     803         952 : RgXY_swap(GEN x, long n, long w)
     804             : {
     805         952 :   GEN z = RgXY_swapspec(x+2, n, w, lgpol(x));
     806         952 :   setvarn(z, varn(x)); return z;
     807             : }
     808             : 
     809             : long
     810         190 : RgXY_degreex(GEN b)
     811             : {
     812         190 :   long deg = 0, i;
     813         190 :   if (!signe(b)) return -1;
     814        1032 :   for (i = 2; i < lg(b); ++i)
     815             :   {
     816         842 :     GEN bi = gel(b, i);
     817         842 :     if (typ(bi) == t_POL)
     818         799 :       deg = maxss(deg, degpol(bi));
     819             :   }
     820         190 :   return deg;
     821             : }
     822             : 
     823             : /* return (x % X^n). Shallow */
     824             : GEN
     825     4562915 : RgXn_red_shallow(GEN a, long n)
     826             : {
     827     4562915 :   long i, L = n+2, l = lg(a);
     828             :   GEN  b;
     829     4562915 :   if (L >= l) return a; /* deg(x) < n */
     830     2972176 :   b = cgetg(L, t_POL); b[1] = a[1];
     831    23208636 :   for (i=2; i<L; i++) gel(b,i) = gel(a,i);
     832     2972176 :   return normalizepol_lg(b,L);
     833             : }
     834             : 
     835             : GEN
     836         357 : RgXnV_red_shallow(GEN x, long n)
     837        1680 : { pari_APPLY_type(t_VEC, RgXn_red_shallow(gel(x,i), n)) }
     838             : 
     839             : /* return (x * X^n). Shallow */
     840             : GEN
     841    75052359 : RgX_shift_shallow(GEN a, long n)
     842             : {
     843    75052359 :   long i, l = lg(a);
     844             :   GEN  b;
     845    75052359 :   if (l == 2 || !n) return a;
     846    55012983 :   l += n;
     847    55012983 :   if (n < 0)
     848             :   {
     849    44790916 :     if (l <= 2) return pol_0(varn(a));
     850    44002732 :     b = cgetg(l, t_POL); b[1] = a[1];
     851    44002732 :     a -= n;
     852   118020451 :     for (i=2; i<l; i++) gel(b,i) = gel(a,i);
     853             :   } else {
     854    10222067 :     b = cgetg(l, t_POL); b[1] = a[1];
     855    10222067 :     a -= n; n += 2;
     856    23305499 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     857   282254342 :     for (   ; i<l; i++) gel(b,i) = gel(a,i);
     858             :   }
     859    54224799 :   return b;
     860             : }
     861             : /* return (x * X^n). */
     862             : GEN
     863     1571386 : RgX_shift(GEN a, long n)
     864             : {
     865     1571386 :   long i, l = lg(a);
     866             :   GEN  b;
     867     1571386 :   if (l == 2 || !n) return RgX_copy(a);
     868     1571106 :   l += n;
     869     1571106 :   if (n < 0)
     870             :   {
     871         630 :     if (l <= 2) return pol_0(varn(a));
     872         588 :     b = cgetg(l, t_POL); b[1] = a[1];
     873         588 :     a -= n;
     874        2275 :     for (i=2; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     875             :   } else {
     876     1570476 :     b = cgetg(l, t_POL); b[1] = a[1];
     877     1570476 :     a -= n; n += 2;
     878     3933163 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     879     4374923 :     for (   ; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     880             :   }
     881     1571064 :   return b;
     882             : }
     883             : 
     884             : GEN
     885      316113 : RgX_rotate_shallow(GEN P, long k, long p)
     886             : {
     887      316113 :   long i, l = lgpol(P);
     888             :   GEN r;
     889      316113 :   if (signe(P)==0)
     890        1365 :     return pol_0(varn(P));
     891      314748 :   r = cgetg(p+2,t_POL); r[1] = P[1];
     892     2090452 :   for(i=0; i<p; i++)
     893             :   {
     894     1775704 :     long s = 2+(i+k)%p;
     895     1775704 :     gel(r,s) = i<l? gel(P,2+i): gen_0;
     896             :   }
     897      314748 :   return RgX_renormalize(r);
     898             : }
     899             : 
     900             : GEN
     901     2975178 : RgX_mulXn(GEN x, long d)
     902             : {
     903             :   pari_sp av;
     904             :   GEN z;
     905             :   long v;
     906     2975178 :   if (d >= 0) return RgX_shift(x, d);
     907     1463667 :   d = -d;
     908     1463667 :   v = RgX_val(x);
     909     1463667 :   if (v >= d) return RgX_shift(x, -d);
     910     1463653 :   av = avma;
     911     1463653 :   z = gred_rfrac_simple(RgX_shift_shallow(x, -v), pol_xn(d - v, varn(x)));
     912     1463653 :   return gerepileupto(av, z);
     913             : }
     914             : 
     915             : long
     916         294 : RgXV_maxdegree(GEN x)
     917             : {
     918         294 :   long d = -1, i, l = lg(x);
     919        2058 :   for (i = 1; i < l; i++)
     920        1764 :     d = maxss(d, degpol(gel(x,i)));
     921         294 :   return d;
     922             : }
     923             : 
     924             : long
     925     2672736 : RgX_val(GEN x)
     926             : {
     927     2672736 :   long i, lx = lg(x);
     928     2672736 :   if (lx == 2) return LONG_MAX;
     929     2699581 :   for (i = 2; i < lx; i++)
     930     2699539 :     if (!isexactzero(gel(x,i))) break;
     931     2672596 :   if (i == lx) return LONG_MAX;/* possible with non-rational zeros */
     932     2672554 :   return i - 2;
     933             : }
     934             : long
     935    39416523 : RgX_valrem(GEN x, GEN *Z)
     936             : {
     937    39416523 :   long v, i, lx = lg(x);
     938    39416523 :   if (lx == 2) { *Z = pol_0(varn(x)); return LONG_MAX; }
     939    74390222 :   for (i = 2; i < lx; i++)
     940    74390219 :     if (!isexactzero(gel(x,i))) break;
     941             :   /* possible with non-rational zeros */
     942    39416526 :   if (i == lx) { *Z = pol_0(varn(x)); return LONG_MAX; }
     943    39416526 :   v = i - 2;
     944    39416526 :   *Z = RgX_shift_shallow(x, -v);
     945    39416523 :   return v;
     946             : }
     947             : long
     948       17794 : RgX_valrem_inexact(GEN x, GEN *Z)
     949             : {
     950             :   long v;
     951       17794 :   if (!signe(x)) { if (Z) *Z = pol_0(varn(x)); return LONG_MAX; }
     952       17787 :   for (v = 0;; v++)
     953       18956 :     if (!gequal0(gel(x,2+v))) break;
     954       17787 :   if (Z) *Z = RgX_shift_shallow(x, -v);
     955       17787 :   return v;
     956             : }
     957             : 
     958             : GEN
     959      109438 : RgXQC_red(GEN x, GEN T)
     960     2051574 : { pari_APPLY_type(t_COL, grem(gel(x,i), T)) }
     961             : 
     962             : GEN
     963        1148 : RgXQV_red(GEN x, GEN T)
     964       27545 : { pari_APPLY_type(t_VEC, grem(gel(x,i), T)) }
     965             : 
     966             : GEN
     967       14826 : RgXQM_red(GEN x, GEN T)
     968      124264 : { pari_APPLY_same(RgXQC_red(gel(x,i), T)) }
     969             : 
     970             : GEN
     971           0 : RgXQM_mul(GEN P, GEN Q, GEN T)
     972             : {
     973           0 :   return RgXQM_red(RgM_mul(P, Q), T);
     974             : }
     975             : 
     976             : GEN
     977      444853 : RgXQX_red(GEN P, GEN T)
     978             : {
     979      444853 :   long i, l = lg(P);
     980      444853 :   GEN Q = cgetg(l, t_POL);
     981      444853 :   Q[1] = P[1];
     982     2976082 :   for (i=2; i<l; i++) gel(Q,i) = grem(gel(P,i), T);
     983      444853 :   return normalizepol_lg(Q, l);
     984             : }
     985             : 
     986             : GEN
     987      235267 : RgX_deriv(GEN x)
     988             : {
     989      235267 :   long i,lx = lg(x)-1;
     990             :   GEN y;
     991             : 
     992      235267 :   if (lx<3) return pol_0(varn(x));
     993      233797 :   y = cgetg(lx,t_POL); gel(y,2) = gcopy(gel(x,3));
     994     1070772 :   for (i=3; i<lx ; i++) gel(y,i) = gmulsg(i-1,gel(x,i+1));
     995      233797 :   y[1] = x[1]; return normalizepol_lg(y,i);
     996             : }
     997             : 
     998             : GEN
     999      496218 : RgX_recipspec_shallow(GEN x, long l, long n)
    1000             : {
    1001             :   long i;
    1002      496218 :   GEN z = cgetg(n+2,t_POL);
    1003      496219 :   z[1] = 0; z += 2;
    1004    15454428 :   for(i=0; i<l; i++)
    1005    14958209 :     gel(z,n-i-1) = gel(x,i);
    1006      624136 :   for(   ; i<n; i++)
    1007      127917 :     gel(z, n-i-1) = gen_0;
    1008      496219 :   return normalizepol_lg(z-2,n+2);
    1009             : }
    1010             : 
    1011             : GEN
    1012       48229 : RgXn_recip_shallow(GEN P, long n)
    1013             : {
    1014       48229 :   GEN Q = RgX_recipspec_shallow(P+2, lgpol(P), n);
    1015       48229 :   setvarn(Q, varn(P));
    1016       48229 :   return Q;
    1017             : }
    1018             : 
    1019             : /* return coefficients s.t x = x_0 X^n + ... + x_n */
    1020             : GEN
    1021        9513 : RgX_recip(GEN x)
    1022             : {
    1023             :   long lx, i, j;
    1024        9513 :   GEN y = cgetg_copy(x, &lx);
    1025       92071 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gcopy(gel(x,j));
    1026        9513 :   return normalizepol_lg(y,lx);
    1027             : }
    1028             : /* shallow version */
    1029             : GEN
    1030      537537 : RgX_recip_shallow(GEN x)
    1031             : {
    1032             :   long lx, i, j;
    1033      537537 :   GEN y = cgetg_copy(x, &lx);
    1034     3729667 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
    1035      537556 :   return y;
    1036             : }
    1037             : /*******************************************************************/
    1038             : /*                                                                 */
    1039             : /*                      ADDITION / SUBTRACTION                     */
    1040             : /*                                                                 */
    1041             : /*******************************************************************/
    1042             : /* same variable */
    1043             : GEN
    1044    34604169 : RgX_add(GEN x, GEN y)
    1045             : {
    1046    34604169 :   long i, lx = lg(x), ly = lg(y);
    1047             :   GEN z;
    1048    34604169 :   if (ly <= lx) {
    1049    30306833 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1050   142051735 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1051    64278859 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1052    30306833 :     z = normalizepol_lg(z, lx);
    1053             :   } else {
    1054     4297336 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1055    29107344 :     for (i=2; i < lx; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1056    15214561 :     for (   ; i < ly; i++) gel(z,i) = gcopy(gel(y,i));
    1057     4297350 :     z = normalizepol_lg(z, ly);
    1058             :   }
    1059    34604186 :   return z;
    1060             : }
    1061             : GEN
    1062    16267859 : RgX_sub(GEN x, GEN y)
    1063             : {
    1064    16267859 :   long i, lx = lg(x), ly = lg(y);
    1065             :   GEN z;
    1066    16267859 :   if (ly <= lx) {
    1067    11103321 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1068    70547957 :     for (i=2; i < ly; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1069    20592537 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1070    11103305 :     z = normalizepol_lg(z, lx);
    1071             :   } else {
    1072     5164538 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1073    22201752 :     for (i=2; i < lx; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1074    13828297 :     for (   ; i < ly; i++) gel(z,i) = gneg(gel(y,i));
    1075     5164537 :     z = normalizepol_lg(z, ly);
    1076             :   }
    1077    16267856 :   return z;
    1078             : }
    1079             : GEN
    1080     1789005 : RgX_neg(GEN x)
    1081             : {
    1082     1789005 :   long i, lx = lg(x);
    1083     1789005 :   GEN y = cgetg(lx, t_POL); y[1] = x[1];
    1084    21830134 :   for (i=2; i<lx; i++) gel(y,i) = gneg(gel(x,i));
    1085     1789005 :   return y;
    1086             : }
    1087             : 
    1088             : GEN
    1089    13144230 : RgX_Rg_add(GEN y, GEN x)
    1090             : {
    1091             :   GEN z;
    1092    13144230 :   long lz = lg(y), i;
    1093    13144230 :   if (lz == 2) return scalarpol(x,varn(y));
    1094    10143131 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1095    10143131 :   gel(z,2) = gadd(gel(y,2),x);
    1096    37729087 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1097             :   /* probably useless unless lz = 3, but cannot be skipped if y is
    1098             :    * an inexact 0 */
    1099    10143131 :   return normalizepol_lg(z,lz);
    1100             : }
    1101             : GEN
    1102        2925 : RgX_Rg_add_shallow(GEN y, GEN x)
    1103             : {
    1104             :   GEN z;
    1105        2925 :   long lz = lg(y), i;
    1106        2925 :   if (lz == 2) return scalarpol(x,varn(y));
    1107        2925 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1108        2925 :   gel(z,2) = gadd(gel(y,2),x);
    1109        5920 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
    1110        2925 :   return z = normalizepol_lg(z,lz);
    1111             : }
    1112             : GEN
    1113       36622 : RgX_Rg_sub(GEN y, GEN x)
    1114             : {
    1115             :   GEN z;
    1116       36622 :   long lz = lg(y), i;
    1117       36622 :   if (lz == 2)
    1118             :   { /* scalarpol(gneg(x),varn(y)) optimized */
    1119        2576 :     long v = varn(y);
    1120        2576 :     if (isrationalzero(x)) return pol_0(v);
    1121          14 :     z = cgetg(3,t_POL);
    1122          14 :     z[1] = gequal0(x)? evalvarn(v)
    1123          14 :                    : evalvarn(v) | evalsigne(1);
    1124          14 :     gel(z,2) = gneg(x); return z;
    1125             :   }
    1126       34046 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1127       34046 :   gel(z,2) = gsub(gel(y,2),x);
    1128      140098 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1129       34046 :   return z = normalizepol_lg(z,lz);
    1130             : }
    1131             : GEN
    1132      636829 : Rg_RgX_sub(GEN x, GEN y)
    1133             : {
    1134             :   GEN z;
    1135      636829 :   long lz = lg(y), i;
    1136      636829 :   if (lz == 2) return scalarpol(x,varn(y));
    1137      635737 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1138      635737 :   gel(z,2) = gsub(x, gel(y,2));
    1139     1294422 :   for(i=3; i<lz; i++) gel(z,i) = gneg(gel(y,i));
    1140      635737 :   return z = normalizepol_lg(z,lz);
    1141             : }
    1142             : /*******************************************************************/
    1143             : /*                                                                 */
    1144             : /*                  KARATSUBA MULTIPLICATION                       */
    1145             : /*                                                                 */
    1146             : /*******************************************************************/
    1147             : #if 0
    1148             : /* to debug Karatsuba-like routines */
    1149             : GEN
    1150             : zx_debug_spec(GEN x, long nx)
    1151             : {
    1152             :   GEN z = cgetg(nx+2,t_POL);
    1153             :   long i;
    1154             :   for (i=0; i<nx; i++) gel(z,i+2) = stoi(x[i]);
    1155             :   z[1] = evalsigne(1); return z;
    1156             : }
    1157             : 
    1158             : GEN
    1159             : RgX_debug_spec(GEN x, long nx)
    1160             : {
    1161             :   GEN z = cgetg(nx+2,t_POL);
    1162             :   long i;
    1163             :   for (i=0; i<nx; i++) z[i+2] = x[i];
    1164             :   z[1] = evalsigne(1); return z;
    1165             : }
    1166             : #endif
    1167             : 
    1168             : /* generic multiplication */
    1169             : GEN
    1170     3451983 : RgX_addspec_shallow(GEN x, GEN y, long nx, long ny)
    1171             : {
    1172             :   GEN z, t;
    1173             :   long i;
    1174     3451983 :   if (nx == ny) {
    1175      726596 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1176     2720065 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1177      726596 :     return normalizepol_lg(z, nx+2);
    1178             :   }
    1179     2725387 :   if (ny < nx) {
    1180     2567658 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1181    11600322 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1182     7934159 :     for (   ; i < nx; i++) gel(t,i) = gel(x,i);
    1183     2567658 :     return normalizepol_lg(z, nx+2);
    1184             :   } else {
    1185      157729 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1186     3120938 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1187      394703 :     for (   ; i < ny; i++) gel(t,i) = gel(y,i);
    1188      157729 :     return normalizepol_lg(z, ny+2);
    1189             :   }
    1190             : }
    1191             : GEN
    1192      208899 : RgX_addspec(GEN x, GEN y, long nx, long ny)
    1193             : {
    1194             :   GEN z, t;
    1195             :   long i;
    1196      208899 :   if (nx == ny) {
    1197       12649 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1198     2143645 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1199       12649 :     return normalizepol_lg(z, nx+2);
    1200             :   }
    1201      196250 :   if (ny < nx) {
    1202      194549 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1203     3214185 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1204     2285423 :     for (   ; i < nx; i++) gel(t,i) = gcopy(gel(x,i));
    1205      194549 :     return normalizepol_lg(z, nx+2);
    1206             :   } else {
    1207        1701 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1208      327040 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1209       12040 :     for (   ; i < ny; i++) gel(t,i) = gcopy(gel(y,i));
    1210        1701 :     return normalizepol_lg(z, ny+2);
    1211             :   }
    1212             : }
    1213             : 
    1214             : /* Return the vector of coefficients of x, where we replace rational 0s by NULL
    1215             :  * [ to speed up basic operation s += x[i]*y[j] ]. We create a proper
    1216             :  * t_VECSMALL, to hold this, which can be left on stack: gerepile
    1217             :  * will not crash on it. The returned vector itself is not a proper GEN,
    1218             :  * we access the coefficients as x[i], i = 0..deg(x) */
    1219             : static GEN
    1220    28880707 : RgXspec_kill0(GEN x, long lx)
    1221             : {
    1222    28880707 :   GEN z = cgetg(lx+1, t_VECSMALL) + 1; /* inhibit gerepile-wise */
    1223             :   long i;
    1224   109670419 :   for (i=0; i <lx; i++)
    1225             :   {
    1226    80789712 :     GEN c = gel(x,i);
    1227    80789712 :     z[i] = (long)(isrationalzero(c)? NULL: c);
    1228             :   }
    1229    28880707 :   return z;
    1230             : }
    1231             : 
    1232             : INLINE GEN
    1233    56429973 : RgX_mulspec_basecase_limb(GEN x, GEN y, long a, long b)
    1234             : {
    1235    56429973 :   pari_sp av = avma;
    1236    56429973 :   GEN s = NULL;
    1237             :   long i;
    1238             : 
    1239   256303011 :   for (i=a; i<b; i++)
    1240   199873038 :     if (gel(y,i) && gel(x,-i))
    1241             :     {
    1242   120213980 :       GEN t = gmul(gel(y,i), gel(x,-i));
    1243   120213980 :       s = s? gadd(s, t): t;
    1244             :     }
    1245    56429973 :   return s? gerepileupto(av, s): gen_0;
    1246             : }
    1247             : 
    1248             : /* assume nx >= ny > 0, return x * y * t^v */
    1249             : static GEN
    1250    10292572 : RgX_mulspec_basecase(GEN x, GEN y, long nx, long ny, long v)
    1251             : {
    1252             :   long i, lz, nz;
    1253             :   GEN z;
    1254             : 
    1255    10292572 :   x = RgXspec_kill0(x,nx);
    1256    10292572 :   y = RgXspec_kill0(y,ny);
    1257    10292572 :   lz = nx + ny + 1; nz = lz-2;
    1258    10292572 :   lz += v;
    1259    10292572 :   z = cgetg(lz, t_POL) + 2; /* x:y:z [i] = term of degree i */
    1260    22066847 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1261    30679483 :   for (i=0; i<ny; i++)gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0, i+1);
    1262    19989227 :   for (  ; i<nx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ny);
    1263    20386911 :   for (  ; i<nz; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-nx+1,ny);
    1264    10292572 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1265             : }
    1266             : 
    1267             : /* return (x * X^d) + y. Assume d > 0 */
    1268             : GEN
    1269     2770094 : RgX_addmulXn_shallow(GEN x0, GEN y0, long d)
    1270             : {
    1271             :   GEN x, y, xd, yd, zd;
    1272             :   long a, lz, nx, ny;
    1273             : 
    1274     2770094 :   if (!signe(x0)) return y0;
    1275     2743064 :   ny = lgpol(y0);
    1276     2743064 :   nx = lgpol(x0);
    1277     2743064 :   zd = (GEN)avma;
    1278     2743064 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1279     2743064 :   if (a <= 0)
    1280             :   {
    1281      228247 :     lz = nx+d+2;
    1282      228247 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1283      892486 :     while (xd > x) gel(--zd,0) = gel(--xd,0);
    1284      228247 :     x = zd + a;
    1285      238178 :     while (zd > x) gel(--zd,0) = gen_0;
    1286             :   }
    1287             :   else
    1288             :   {
    1289     2514817 :     xd = new_chunk(d); yd = y+d;
    1290     2514817 :     x = RgX_addspec_shallow(x,yd, nx,a);
    1291     2514817 :     lz = (a>nx)? ny+2: lg(x)+d;
    1292    17467706 :     x += 2; while (xd > x) *--zd = *--xd;
    1293             :   }
    1294     8582284 :   while (yd > y) *--zd = *--yd;
    1295     2743064 :   *--zd = x0[1];
    1296     2743064 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1297             : }
    1298             : GEN
    1299      501171 : RgX_addmulXn(GEN x0, GEN y0, long d)
    1300             : {
    1301             :   GEN x, y, xd, yd, zd;
    1302             :   long a, lz, nx, ny;
    1303             : 
    1304      501171 :   if (!signe(x0)) return RgX_copy(y0);
    1305      500401 :   nx = lgpol(x0);
    1306      500401 :   ny = lgpol(y0);
    1307      500401 :   zd = (GEN)avma;
    1308      500401 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1309      500401 :   if (a <= 0)
    1310             :   {
    1311      291502 :     lz = nx+d+2;
    1312      291502 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1313     4259859 :     while (xd > x) gel(--zd,0) = gcopy(gel(--xd,0));
    1314      291502 :     x = zd + a;
    1315      757381 :     while (zd > x) gel(--zd,0) = gen_0;
    1316             :   }
    1317             :   else
    1318             :   {
    1319      208899 :     xd = new_chunk(d); yd = y+d;
    1320      208899 :     x = RgX_addspec(x,yd, nx,a);
    1321      208899 :     lz = (a>nx)? ny+2: lg(x)+d;
    1322     7786083 :     x += 2; while (xd > x) *--zd = *--xd;
    1323             :   }
    1324     2542887 :   while (yd > y) gel(--zd,0) = gcopy(gel(--yd,0));
    1325      500401 :   *--zd = x0[1];
    1326      500401 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1327             : }
    1328             : 
    1329             : /* return x * y mod t^n */
    1330             : static GEN
    1331     4143023 : RgXn_mul_basecase(GEN x, GEN y, long n)
    1332             : {
    1333     4143023 :   long i, lz = n+2, lx = lgpol(x), ly = lgpol(y);
    1334             :   GEN z;
    1335     4143023 :   if (lx < 0) return pol_0(varn(x));
    1336     4143023 :   if (ly < 0) return pol_0(varn(x));
    1337     4143023 :   z = cgetg(lz, t_POL) + 2;
    1338     4143023 :   x+=2; if (lx > n) lx = n;
    1339     4143023 :   y+=2; if (ly > n) ly = n;
    1340     4143023 :   z[-1] = x[-1];
    1341     4143023 :   if (ly > lx) { swap(x,y); lswap(lx,ly); }
    1342     4143023 :   x = RgXspec_kill0(x, lx);
    1343     4143023 :   y = RgXspec_kill0(y, ly);
    1344             :   /* x:y:z [i] = term of degree i */
    1345    18228515 :   for (i=0;i<ly; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,i+1);
    1346     6261720 :   for (  ; i<lx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ly);
    1347     4190902 :   for (  ; i<n; i++)  gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-lx+1,ly);
    1348     4143023 :   return normalizepol_lg(z - 2, lz);
    1349             : }
    1350             : /* Mulders / Karatsuba product f*g mod t^n (Hanrot-Zimmermann variant) */
    1351             : static GEN
    1352     5152273 : RgXn_mul2(GEN f, GEN g, long n)
    1353             : {
    1354     5152273 :   pari_sp av = avma;
    1355             :   GEN fe,fo, ge,go, l,h,m;
    1356             :   long n0, n1;
    1357     5152273 :   if (degpol(f) + degpol(g) < n) return RgX_mul(f,g);
    1358     4174425 :   if (n < 80) return RgXn_mul_basecase(f,g,n);
    1359       31402 :   n0 = n>>1; n1 = n-n0;
    1360       31402 :   RgX_even_odd(f, &fe, &fo);
    1361       31402 :   RgX_even_odd(g, &ge, &go);
    1362       31402 :   l = RgXn_mul(fe,ge,n1);
    1363       31402 :   h = RgXn_mul(fo,go,n0);
    1364       31402 :   m = RgX_sub(RgXn_mul(RgX_add(fe,fo),RgX_add(ge,go),n0), RgX_add(l,h));
    1365             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1366             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1367       31402 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1368             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1369       31402 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1370       31402 :   m = RgX_inflate(m,2);
    1371             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1372       31402 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1373       31402 :   h = RgX_inflate(h,2);
    1374       31402 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1375       31402 :   return gerepileupto(av, h);
    1376             : }
    1377             : /* (f*g) \/ x^n */
    1378             : static GEN
    1379      885463 : RgX_mulhigh_i2(GEN f, GEN g, long n)
    1380             : {
    1381      885463 :   long d = degpol(f)+degpol(g) + 1 - n;
    1382             :   GEN h;
    1383      885463 :   if (d <= 2) return RgX_shift_shallow(RgX_mul(f,g), -n);
    1384       21310 :   h = RgX_recip_shallow(RgXn_mul(RgX_recip_shallow(f),
    1385             :                                  RgX_recip_shallow(g), d));
    1386       21310 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1387             : }
    1388             : 
    1389             : /* (f*g) \/ x^n */
    1390             : static GEN
    1391           0 : RgX_sqrhigh_i2(GEN f, long n)
    1392             : {
    1393           0 :   long d = 2*degpol(f)+ 1 - n;
    1394             :   GEN h;
    1395           0 :   if (d <= 2) return RgX_shift_shallow(RgX_sqr(f), -n);
    1396           0 :   h = RgX_recip_shallow(RgXn_sqr(RgX_recip_shallow(f), d));
    1397           0 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1398             : }
    1399             : 
    1400             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
    1401             :  * b+2 were sent instead. na, nb = number of terms of a, b.
    1402             :  * Only c, c0, c1, c2 are genuine GEN.
    1403             :  */
    1404             : GEN
    1405    11003729 : RgX_mulspec(GEN a, GEN b, long na, long nb)
    1406             : {
    1407             :   GEN a0, c, c0;
    1408    11003729 :   long n0, n0a, i, v = 0;
    1409             :   pari_sp av;
    1410             : 
    1411    17409317 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v++; }
    1412    16907082 :   while (nb && isrationalzero(gel(b,0))) { b++; nb--; v++; }
    1413    11003729 :   if (na < nb) swapspec(a,b, na,nb);
    1414    11003729 :   if (!nb) return pol_0(0);
    1415             : 
    1416    10762302 :   if (nb < RgX_MUL_LIMIT) return RgX_mulspec_basecase(a,b,na,nb, v);
    1417      469730 :   RgX_shift_inplace_init(v);
    1418      469730 :   i = (na>>1); n0 = na-i; na = i;
    1419      469730 :   av = avma; a0 = a+n0; n0a = n0;
    1420     1325230 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1421             : 
    1422      469730 :   if (nb > n0)
    1423             :   {
    1424             :     GEN b0,c1,c2;
    1425             :     long n0b;
    1426             : 
    1427      468583 :     nb -= n0; b0 = b+n0; n0b = n0;
    1428     1436546 :     while (n0b && isrationalzero(gel(b,n0b-1))) n0b--;
    1429      468583 :     c = RgX_mulspec(a,b,n0a,n0b);
    1430      468583 :     c0 = RgX_mulspec(a0,b0, na,nb);
    1431             : 
    1432      468583 :     c2 = RgX_addspec_shallow(a0,a, na,n0a);
    1433      468583 :     c1 = RgX_addspec_shallow(b0,b, nb,n0b);
    1434             : 
    1435      468583 :     c1 = RgX_mulspec(c1+2,c2+2, lgpol(c1),lgpol(c2));
    1436      468583 :     c2 = RgX_sub(c1, RgX_add(c0,c));
    1437      468583 :     c0 = RgX_addmulXn_shallow(c0, c2, n0);
    1438             :   }
    1439             :   else
    1440             :   {
    1441        1147 :     c = RgX_mulspec(a,b,n0a,nb);
    1442        1147 :     c0 = RgX_mulspec(a0,b,na,nb);
    1443             :   }
    1444      469730 :   c0 = RgX_addmulXn(c0,c,n0);
    1445      469730 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1446             : }
    1447             : 
    1448             : INLINE GEN
    1449       49591 : RgX_sqrspec_basecase_limb(GEN x, long a, long i)
    1450             : {
    1451       49591 :   pari_sp av = avma;
    1452       49591 :   GEN s = NULL;
    1453       49591 :   long j, l = (i+1)>>1;
    1454      127995 :   for (j=a; j<l; j++)
    1455             :   {
    1456       78404 :     GEN xj = gel(x,j), xx = gel(x,i-j);
    1457       78404 :     if (xj && xx)
    1458             :     {
    1459       70662 :       GEN t = gmul(xj, xx);
    1460       70662 :       s = s? gadd(s, t): t;
    1461             :     }
    1462             :   }
    1463       49591 :   if (s) s = gshift(s,1);
    1464       49591 :   if ((i&1) == 0)
    1465             :   {
    1466       29554 :     GEN t = gel(x, i>>1);
    1467       29554 :     if (t) {
    1468       27118 :       t = gsqr(t);
    1469       27118 :       s = s? gadd(s, t): t;
    1470             :     }
    1471             :   }
    1472       49591 :   return s? gerepileupto(av,s): gen_0;
    1473             : }
    1474             : static GEN
    1475        9517 : RgX_sqrspec_basecase(GEN x, long nx, long v)
    1476             : {
    1477             :   long i, lz, nz;
    1478             :   GEN z;
    1479             : 
    1480        9517 :   if (!nx) return pol_0(0);
    1481        9517 :   x = RgXspec_kill0(x,nx);
    1482        9517 :   lz = (nx << 1) + 1, nz = lz-2;
    1483        9517 :   lz += v;
    1484        9517 :   z = cgetg(lz,t_POL) + 2;
    1485       18197 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1486       39071 :   for (i=0; i<nx; i++)gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1487       29554 :   for (  ; i<nz; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, i-nx+1, i);
    1488        9517 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1489             : }
    1490             : /* return x^2 mod t^n */
    1491             : static GEN
    1492           0 : RgXn_sqr_basecase(GEN x, long n)
    1493             : {
    1494           0 :   long i, lz = n+2, lx = lgpol(x);
    1495             :   GEN z;
    1496           0 :   if (lx < 0) return pol_0(varn(x));
    1497           0 :   z = cgetg(lz, t_POL);
    1498           0 :   z[1] = x[1];
    1499           0 :   x+=2; if (lx > n) lx = n;
    1500           0 :   x = RgXspec_kill0(x,lx);
    1501           0 :   z+=2;/* x:z [i] = term of degree i */
    1502           0 :   for (i=0;i<lx; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1503           0 :   for (  ; i<n; i++)  gel(z,i) = RgX_sqrspec_basecase_limb(x, i-lx+1, i);
    1504           0 :   z -= 2; return normalizepol_lg(z, lz);
    1505             : }
    1506             : /* Mulders / Karatsuba product f^2 mod t^n (Hanrot-Zimmermann variant) */
    1507             : static GEN
    1508         245 : RgXn_sqr2(GEN f, long n)
    1509             : {
    1510         245 :   pari_sp av = avma;
    1511             :   GEN fe,fo, l,h,m;
    1512             :   long n0, n1;
    1513         245 :   if (2*degpol(f) < n) return RgX_sqr_i(f);
    1514           0 :   if (n < 80) return RgXn_sqr_basecase(f,n);
    1515           0 :   n0 = n>>1; n1 = n-n0;
    1516           0 :   RgX_even_odd(f, &fe, &fo);
    1517           0 :   l = RgXn_sqr(fe,n1);
    1518           0 :   h = RgXn_sqr(fo,n0);
    1519           0 :   m = RgX_sub(RgXn_sqr(RgX_add(fe,fo),n0), RgX_add(l,h));
    1520             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1521             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1522           0 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1523             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1524           0 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1525           0 :   m = RgX_inflate(m,2);
    1526             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1527           0 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1528           0 :   h = RgX_inflate(h,2);
    1529           0 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1530           0 :   return gerepileupto(av, h);
    1531             : }
    1532             : GEN
    1533        9556 : RgX_sqrspec(GEN a, long na)
    1534             : {
    1535             :   GEN a0, c, c0, c1;
    1536        9556 :   long n0, n0a, i, v = 0;
    1537             :   pari_sp av;
    1538             : 
    1539       13896 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v += 2; }
    1540        9556 :   if (na<RgX_SQR_LIMIT) return RgX_sqrspec_basecase(a, na, v);
    1541          39 :   RgX_shift_inplace_init(v);
    1542          39 :   i = (na>>1); n0 = na-i; na = i;
    1543          39 :   av = avma; a0 = a+n0; n0a = n0;
    1544          39 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1545             : 
    1546          39 :   c = RgX_sqrspec(a,n0a);
    1547          39 :   c0 = RgX_sqrspec(a0,na);
    1548          39 :   c1 = gmul2n(RgX_mulspec(a0,a, na,n0a), 1);
    1549          39 :   c0 = RgX_addmulXn_shallow(c0,c1, n0);
    1550          39 :   c0 = RgX_addmulXn(c0,c,n0);
    1551          39 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1552             : }
    1553             : 
    1554             : /* (X^a + A)(X^b + B) - X^(a+b), where deg A < a, deg B < b */
    1555             : GEN
    1556      783886 : RgX_mul_normalized(GEN A, long a, GEN B, long b)
    1557             : {
    1558      783886 :   GEN z = RgX_mul(A, B);
    1559      783886 :   if (a < b)
    1560        5446 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(A, B, b-a), z, a);
    1561      778440 :   else if (a > b)
    1562      499227 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(B, A, a-b), z, b);
    1563             :   else
    1564      279213 :     z = RgX_addmulXn_shallow(RgX_add(A, B), z, a);
    1565      783886 :   return z;
    1566             : }
    1567             : 
    1568             : GEN
    1569     9595647 : RgX_mul_i(GEN x, GEN y)
    1570             : {
    1571     9595647 :   GEN z = RgX_mulspec(x+2, y+2, lgpol(x), lgpol(y));
    1572     9595647 :   setvarn(z, varn(x)); return z;
    1573             : }
    1574             : 
    1575             : GEN
    1576        9478 : RgX_sqr_i(GEN x)
    1577             : {
    1578        9478 :   GEN z = RgX_sqrspec(x+2, lgpol(x));
    1579        9478 :   setvarn(z,varn(x)); return z;
    1580             : }
    1581             : 
    1582             : /*******************************************************************/
    1583             : /*                                                                 */
    1584             : /*                               DIVISION                          */
    1585             : /*                                                                 */
    1586             : /*******************************************************************/
    1587             : GEN
    1588      674869 : RgX_Rg_divexact(GEN x, GEN y) {
    1589      674869 :   long i, lx = lg(x);
    1590             :   GEN z;
    1591      674869 :   if (lx == 2) return gcopy(x);
    1592      639881 :   switch(typ(y))
    1593             :   {
    1594      620723 :     case t_INT:
    1595      620723 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1596      587054 :       break;
    1597        3955 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1598             :   }
    1599      602257 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1600     5951477 :   for (i=2; i<lx; i++) gel(z,i) = gdivexact(gel(x,i),y);
    1601      602257 :   return z;
    1602             : }
    1603             : GEN
    1604    24716250 : RgX_Rg_div(GEN x, GEN y) {
    1605    24716250 :   long i, lx = lg(x);
    1606             :   GEN z;
    1607    24716250 :   if (lx == 2) return gcopy(x);
    1608    24497363 :   switch(typ(y))
    1609             :   {
    1610    15372622 :     case t_INT:
    1611    15372622 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1612     1133511 :       break;
    1613        4186 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1614             :   }
    1615    10254066 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1616    37897074 :   for (i=2; i<lx; i++) gel(z,i) = gdiv(gel(x,i),y);
    1617    10254066 :   return normalizepol_lg(z, lx);
    1618             : }
    1619             : GEN
    1620       13377 : RgX_normalize(GEN x)
    1621             : {
    1622       13377 :   GEN z, d = NULL;
    1623       13377 :   long i, n = lg(x)-1;
    1624       13377 :   for (i = n; i > 1; i--) { d = gel(x,i); if (!gequal0(d)) break; }
    1625       13377 :   if (i == 1) return pol_0(varn(x));
    1626       13377 :   if (i == n && isint1(d)) return x;
    1627        5159 :   n = i; z = cgetg(n+1, t_POL); z[1] = x[1];
    1628       11802 :   for (i=2; i<n; i++) gel(z,i) = gdiv(gel(x,i),d);
    1629        5159 :   gel(z,n) = Rg_get_1(d); return z;
    1630             : }
    1631             : GEN
    1632        2625 : RgX_divs(GEN x, long y) {
    1633             :   long i, lx;
    1634        2625 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1635       10318 :   for (i=2; i<lx; i++) gel(z,i) = gdivgs(gel(x,i),y);
    1636        2625 :   return normalizepol_lg(z, lx);
    1637             : }
    1638             : GEN
    1639       60189 : RgX_div_by_X_x(GEN a, GEN x, GEN *r)
    1640             : {
    1641       60189 :   long l = lg(a), i;
    1642             :   GEN a0, z0, z;
    1643             : 
    1644       60189 :   if (l <= 3)
    1645             :   {
    1646           0 :     if (r) *r = l == 2? gen_0: gcopy(gel(a,2));
    1647           0 :     return pol_0(0);
    1648             :   }
    1649       60189 :   z = cgetg(l-1, t_POL);
    1650       60189 :   z[1] = a[1];
    1651       60189 :   a0 = a + l-1;
    1652       60189 :   z0 = z + l-2; *z0 = *a0--;
    1653      954726 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1654             :   {
    1655      894537 :     GEN t = gadd(gel(a0--,0), gmul(x, gel(z0--,0)));
    1656      894537 :     gel(z0,0) = t;
    1657             :   }
    1658       60189 :   if (r) *r = gadd(gel(a0,0), gmul(x, gel(z0,0)));
    1659       60189 :   return z;
    1660             : }
    1661             : /* Polynomial division x / y:
    1662             :  *   if pr = ONLY_REM return remainder, otherwise return quotient
    1663             :  *   if pr = ONLY_DIVIDES return quotient if division is exact, else NULL
    1664             :  *   if pr != NULL set *pr to remainder, as the last object on stack */
    1665             : /* assume, typ(x) = typ(y) = t_POL, same variable */
    1666             : static GEN
    1667    12858869 : RgX_divrem_i(GEN x, GEN y, GEN *pr)
    1668             : {
    1669             :   pari_sp avy, av, av1;
    1670             :   long dx,dy,dz,i,j,sx,lr;
    1671             :   GEN z,p1,p2,rem,y_lead,mod,p;
    1672             :   GEN (*f)(GEN,GEN);
    1673             : 
    1674    12858869 :   if (!signe(y)) pari_err_INV("RgX_divrem",y);
    1675             : 
    1676    12858869 :   dy = degpol(y);
    1677    12858869 :   y_lead = gel(y,dy+2);
    1678    12858869 :   if (gequal0(y_lead)) /* normalize denominator if leading term is 0 */
    1679             :   {
    1680           0 :     pari_warn(warner,"normalizing a polynomial with 0 leading term");
    1681           0 :     for (dy--; dy>=0; dy--)
    1682             :     {
    1683           0 :       y_lead = gel(y,dy+2);
    1684           0 :       if (!gequal0(y_lead)) break;
    1685             :     }
    1686             :   }
    1687    12858869 :   if (!dy) /* y is constant */
    1688             :   {
    1689        3051 :     if (pr == ONLY_REM) return pol_0(varn(x));
    1690        3051 :     z = RgX_Rg_div(x, y_lead);
    1691        3051 :     if (pr == ONLY_DIVIDES) return z;
    1692        1168 :     if (pr) *pr = pol_0(varn(x));
    1693        1168 :     return z;
    1694             :   }
    1695    12855818 :   dx = degpol(x);
    1696    12855818 :   if (dx < dy)
    1697             :   {
    1698      710723 :     if (pr == ONLY_REM) return RgX_copy(x);
    1699      371761 :     if (pr == ONLY_DIVIDES) return signe(x)? NULL: pol_0(varn(x));
    1700      371740 :     z = pol_0(varn(x));
    1701      371740 :     if (pr) *pr = RgX_copy(x);
    1702      371740 :     return z;
    1703             :   }
    1704             : 
    1705             :   /* x,y in R[X], y non constant */
    1706    12145095 :   av = avma;
    1707    12145095 :   p = NULL;
    1708    12145095 :   if (RgX_is_FpX(x, &p) && RgX_is_FpX(y, &p) && p)
    1709             :   {
    1710      225841 :     z = FpX_divrem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p, pr);
    1711      225841 :     if (!z) return gc_NULL(av);
    1712      225841 :     z = FpX_to_mod(z, p);
    1713      225841 :     if (!pr || pr == ONLY_REM || pr == ONLY_DIVIDES)
    1714      121681 :       return gerepileupto(av, z);
    1715      104160 :     *pr = FpX_to_mod(*pr, p);
    1716      104160 :     gerepileall(av, 2, pr, &z);
    1717      104160 :     return z;
    1718             :   }
    1719    11919254 :   switch(typ(y_lead))
    1720             :   {
    1721           0 :     case t_REAL:
    1722           0 :       y_lead = ginv(y_lead);
    1723           0 :       f = gmul; mod = NULL;
    1724           0 :       break;
    1725        1957 :     case t_INTMOD:
    1726        1957 :     case t_POLMOD: y_lead = ginv(y_lead);
    1727        1957 :       f = gmul; mod = gmodulo(gen_1, gel(y_lead,1));
    1728        1957 :       break;
    1729    11917297 :     default: if (gequal1(y_lead)) y_lead = NULL;
    1730    11917297 :       f = gdiv; mod = NULL;
    1731             :   }
    1732             : 
    1733    11919254 :   if (y_lead == NULL)
    1734    10114172 :     p2 = gel(x,dx+2);
    1735             :   else {
    1736             :     for(;;) {
    1737     1805082 :       p2 = f(gel(x,dx+2),y_lead);
    1738     1805082 :       p2 = simplify_shallow(p2);
    1739     1805082 :       if (!isexactzero(p2) || (--dx < 0)) break;
    1740             :     }
    1741     1805082 :     if (dx < dy) /* leading coeff of x was in fact zero */
    1742             :     {
    1743           0 :       if (pr == ONLY_DIVIDES) {
    1744           0 :         set_avma(av);
    1745           0 :         return (dx < 0)? pol_0(varn(x)) : NULL;
    1746             :       }
    1747           0 :       if (pr == ONLY_REM)
    1748             :       {
    1749           0 :         if (dx < 0)
    1750           0 :           return gerepilecopy(av, scalarpol(p2, varn(x)));
    1751             :         else
    1752             :         {
    1753             :           GEN t;
    1754           0 :           set_avma(av);
    1755           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1756           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1757           0 :           return t;
    1758             :         }
    1759             :       }
    1760           0 :       if (pr) /* cf ONLY_REM above */
    1761             :       {
    1762           0 :         if (dx < 0)
    1763             :         {
    1764           0 :           p2 = gclone(p2);
    1765           0 :           set_avma(av);
    1766           0 :           z = pol_0(varn(x));
    1767           0 :           x = scalarpol(p2, varn(x));
    1768           0 :           gunclone(p2);
    1769             :         }
    1770             :         else
    1771             :         {
    1772             :           GEN t;
    1773           0 :           set_avma(av);
    1774           0 :           z = pol_0(varn(x));
    1775           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1776           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1777           0 :           x = t;
    1778             :         }
    1779           0 :         *pr = x;
    1780             :       }
    1781             :       else
    1782             :       {
    1783           0 :         set_avma(av);
    1784           0 :         z = pol_0(varn(x));
    1785             :       }
    1786           0 :       return z;
    1787             :     }
    1788             :   }
    1789             :   /* dx >= dy */
    1790    11919254 :   avy = avma;
    1791    11919254 :   dz = dx-dy;
    1792    11919254 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1793    11919254 :   x += 2;
    1794    11919254 :   z += 2;
    1795    11919254 :   y += 2;
    1796    11919254 :   gel(z,dz) = gcopy(p2);
    1797             : 
    1798    31111072 :   for (i=dx-1; i>=dy; i--)
    1799             :   {
    1800    19191818 :     av1=avma; p1=gel(x,i);
    1801  1069525510 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1802    19191818 :     if (y_lead) p1 = simplify(f(p1,y_lead));
    1803             : 
    1804    19191818 :     if (isrationalzero(p1)) { set_avma(av1); p1 = gen_0; }
    1805             :     else
    1806     9641779 :       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
    1807    19191818 :     gel(z,i-dy) = p1;
    1808             :   }
    1809    11919254 :   if (!pr) return gerepileupto(av,z-2);
    1810             : 
    1811     5231781 :   rem = (GEN)avma; av1 = (pari_sp)new_chunk(dx+3);
    1812     5231781 :   for (sx=0; ; i--)
    1813             :   {
    1814      348887 :     p1 = gel(x,i);
    1815             :     /* we always enter this loop at least once */
    1816    15407551 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1817     5580668 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1818     5580668 :     if (!gequal0(p1)) { sx = 1; break; } /* remainder is non-zero */
    1819     3421663 :     if (!isexactzero(p1)) break;
    1820     3404914 :     if (!i) break;
    1821      348887 :     set_avma(av1);
    1822             :   }
    1823     5231781 :   if (pr == ONLY_DIVIDES)
    1824             :   {
    1825        1659 :     if (sx) return gc_NULL(av);
    1826        1652 :     set_avma((pari_sp)rem); return gerepileupto(av,z-2);
    1827             :   }
    1828     5230122 :   lr=i+3; rem -= lr;
    1829     5230122 :   if (avma==av1) { set_avma((pari_sp)rem); p1 = gcopy(p1); }
    1830     5186552 :   else p1 = gerepileupto((pari_sp)rem,p1);
    1831     5230122 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1832     5230122 :   rem[1] = z[-1];
    1833     5230122 :   rem += 2;
    1834     5230122 :   gel(rem,i) = p1;
    1835     6828672 :   for (i--; i>=0; i--)
    1836             :   {
    1837     1598550 :     av1=avma; p1 = gel(x,i);
    1838     5072768 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1839     1598550 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1840     1598550 :     gel(rem,i) = avma==av1? gcopy(p1):gerepileupto(av1,p1);
    1841             :   }
    1842     5230122 :   rem -= 2;
    1843     5230122 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1844     5230122 :   if (pr == ONLY_REM) return gerepileupto(av,rem);
    1845     4167031 :   z -= 2;
    1846             :   {
    1847     4167031 :     GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
    1848     4167031 :     gerepilemanysp(av,avy,gptr,2); *pr = rem; return z;
    1849             :   }
    1850             : }
    1851             : 
    1852             : GEN
    1853    11456816 : RgX_divrem(GEN x, GEN y, GEN *pr)
    1854             : {
    1855    11456816 :   if (pr == ONLY_REM) return RgX_rem(x, y);
    1856    11456816 :   return RgX_divrem_i(x, y, pr);
    1857             : }
    1858             : 
    1859             : /* x and y in (R[Y]/T)[X]  (lifted), T in R[Y]. y preferably monic */
    1860             : GEN
    1861      110929 : RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)
    1862             : {
    1863             :   long vx, dx, dy, dz, i, j, sx, lr;
    1864             :   pari_sp av0, av, tetpil;
    1865             :   GEN z,p1,rem,lead;
    1866             : 
    1867      110929 :   if (!signe(y)) pari_err_INV("RgXQX_divrem",y);
    1868      110929 :   vx = varn(x);
    1869      110929 :   dx = degpol(x);
    1870      110929 :   dy = degpol(y);
    1871      110929 :   if (dx < dy)
    1872             :   {
    1873       25760 :     if (pr)
    1874             :     {
    1875       25760 :       av0 = avma; x = RgXQX_red(x, T);
    1876       25760 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: gen_0; }
    1877       25760 :       if (pr == ONLY_REM) return x;
    1878           0 :       *pr = x;
    1879             :     }
    1880           0 :     return pol_0(vx);
    1881             :   }
    1882       85169 :   lead = leading_coeff(y);
    1883       85169 :   if (!dy) /* y is constant */
    1884             :   {
    1885         546 :     if (pr && pr != ONLY_DIVIDES)
    1886             :     {
    1887           0 :       if (pr == ONLY_REM) return pol_0(vx);
    1888           0 :       *pr = pol_0(vx);
    1889             :     }
    1890         546 :     if (gequal1(lead)) return RgX_copy(x);
    1891           0 :     av0 = avma; x = gmul(x, ginvmod(lead,T)); tetpil = avma;
    1892           0 :     return gerepile(av0,tetpil,RgXQX_red(x,T));
    1893             :   }
    1894       84623 :   av0 = avma; dz = dx-dy;
    1895       84623 :   lead = gequal1(lead)? NULL: gclone(ginvmod(lead,T));
    1896       84623 :   set_avma(av0);
    1897       84623 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1898       84623 :   x += 2; y += 2; z += 2;
    1899             : 
    1900       84623 :   p1 = gel(x,dx); av = avma;
    1901       84623 :   gel(z,dz) = lead? gerepileupto(av, grem(gmul(p1,lead), T)): gcopy(p1);
    1902      353472 :   for (i=dx-1; i>=dy; i--)
    1903             :   {
    1904      268849 :     av=avma; p1=gel(x,i);
    1905     1563828 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1906      268849 :     if (lead) p1 = gmul(grem(p1, T), lead);
    1907      268849 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil, grem(p1, T));
    1908             :   }
    1909       84623 :   if (!pr) { guncloneNULL(lead); return z-2; }
    1910             : 
    1911       83377 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    1912       83377 :   for (sx=0; ; i--)
    1913             :   {
    1914       33745 :     p1 = gel(x,i);
    1915      409040 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1916      117121 :     tetpil=avma; p1 = grem(p1, T); if (!gequal0(p1)) { sx = 1; break; }
    1917       42866 :     if (!i) break;
    1918       33745 :     set_avma(av);
    1919             :   }
    1920       83377 :   if (pr == ONLY_DIVIDES)
    1921             :   {
    1922        6167 :     guncloneNULL(lead);
    1923        6167 :     if (sx) return gc_NULL(av0);
    1924        6027 :     return gc_const((pari_sp)rem, z-2);
    1925             :   }
    1926       77210 :   lr=i+3; rem -= lr;
    1927       77210 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1928       77210 :   rem[1] = z[-1];
    1929       77210 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    1930       77210 :   rem += 2; gel(rem,i) = p1;
    1931      164914 :   for (i--; i>=0; i--)
    1932             :   {
    1933       87704 :     av=avma; p1 = gel(x,i);
    1934      299693 :     for (j=0; j<=i && j<=dz; j++)
    1935      211989 :       p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1936       87704 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, grem(p1, T));
    1937             :   }
    1938       77210 :   rem -= 2;
    1939       77210 :   guncloneNULL(lead);
    1940       77210 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1941       77210 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    1942         126 :   *pr = rem; return z-2;
    1943             : }
    1944             : 
    1945             : /*******************************************************************/
    1946             : /*                                                                 */
    1947             : /*                        PSEUDO-DIVISION                          */
    1948             : /*                                                                 */
    1949             : /*******************************************************************/
    1950             : INLINE GEN
    1951      504857 : rem(GEN c, GEN T)
    1952             : {
    1953      504857 :   if (T && typ(c) == t_POL && varn(c) == varn(T)) c = RgX_rem(c, T);
    1954      504857 :   return c;
    1955             : }
    1956             : 
    1957             : /* x, y, are ZYX, lc(y) is an integer, T is a ZY */
    1958             : int
    1959        8486 : ZXQX_dvd(GEN x, GEN y, GEN T)
    1960             : {
    1961             :   long dx, dy, dz, i, p, T_ismonic;
    1962        8486 :   pari_sp av = avma, av2;
    1963             :   GEN y_lead;
    1964             : 
    1965        8486 :   if (!signe(y)) pari_err_INV("ZXQX_dvd",y);
    1966        8486 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1967        8486 :   if (typ(y_lead) == t_POL) y_lead = gel(y_lead, 2); /* t_INT */
    1968             :   /* if monic, no point in using pseudo-division */
    1969        8486 :   if (gequal1(y_lead)) return signe(RgXQX_rem(x, y, T)) == 0;
    1970        5392 :   T_ismonic = gequal1(leading_coeff(T));
    1971        5392 :   dx = degpol(x);
    1972        5392 :   if (dx < dy) return !signe(x);
    1973        5392 :   (void)new_chunk(2);
    1974        5392 :   x = RgX_recip_shallow(x)+2;
    1975        5392 :   y = RgX_recip_shallow(y)+2;
    1976             :   /* pay attention to sparse divisors */
    1977       11280 :   for (i = 1; i <= dy; i++)
    1978        5888 :     if (!signe(gel(y,i))) gel(y,i) = NULL;
    1979        5392 :   dz = dx-dy; p = dz+1;
    1980        5392 :   av2 = avma;
    1981             :   for (;;)
    1982       34299 :   {
    1983       39691 :     GEN m, x0 = gel(x,0), y0 = y_lead, cx = content(x0);
    1984       39691 :     x0 = gneg(x0); p--;
    1985       39691 :     m = gcdii(cx, y0);
    1986       39691 :     if (!equali1(m))
    1987             :     {
    1988       37834 :       x0 = gdiv(x0, m);
    1989       37834 :       y0 = diviiexact(y0, m);
    1990       37834 :       if (equali1(y0)) y0 = NULL;
    1991             :     }
    1992       83095 :     for (i=1; i<=dy; i++)
    1993             :     {
    1994       43404 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1995       43404 :       if (gel(y,i)) c = gadd(c, gmul(x0,gel(y,i)));
    1996       43404 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1997       43404 :       gel(x,i) = c;
    1998             :     }
    1999      422575 :     for (   ; i<=dx; i++)
    2000             :     {
    2001      382884 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    2002      382884 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    2003      382884 :       gel(x,i) = c;
    2004             :     }
    2005       45487 :     do { x++; dx--; } while (dx >= 0 && !signe(gel(x,0)));
    2006       39691 :     if (dx < dy) break;
    2007       34299 :     if (gc_needed(av2,1))
    2008             :     {
    2009          28 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZXQX_dvd dx = %ld >= %ld",dx,dy);
    2010          28 :       gerepilecoeffs(av2,x,dx+1);
    2011             :     }
    2012             :   }
    2013        5392 :   return gc_bool(av, dx < 0);
    2014             : }
    2015             : 
    2016             : /* T either NULL or a t_POL. */
    2017             : GEN
    2018       28122 : RgXQX_pseudorem(GEN x, GEN y, GEN T)
    2019             : {
    2020       28122 :   long vx = varn(x), dx, dy, dz, i, lx, p;
    2021       28122 :   pari_sp av = avma, av2;
    2022             :   GEN y_lead;
    2023             : 
    2024       28122 :   if (!signe(y)) pari_err_INV("RgXQX_pseudorem",y);
    2025       28122 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2026             :   /* if monic, no point in using pseudo-division */
    2027       28122 :   if (gequal1(y_lead)) return T? RgXQX_rem(x, y, T): RgX_rem(x, y);
    2028       24951 :   dx = degpol(x);
    2029       24951 :   if (dx < dy) return RgX_copy(x);
    2030       24944 :   (void)new_chunk(2);
    2031       24944 :   x = RgX_recip_shallow(x)+2;
    2032       24944 :   y = RgX_recip_shallow(y)+2;
    2033             :   /* pay attention to sparse divisors */
    2034       79315 :   for (i = 1; i <= dy; i++)
    2035       54371 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2036       24944 :   dz = dx-dy; p = dz+1;
    2037       24944 :   av2 = avma;
    2038             :   for (;;)
    2039             :   {
    2040       90647 :     gel(x,0) = gneg(gel(x,0)); p--;
    2041      297431 :     for (i=1; i<=dy; i++)
    2042             :     {
    2043      206784 :       GEN c = gmul(y_lead, gel(x,i));
    2044      206784 :       if (gel(y,i)) c = gadd(c, gmul(gel(x,0),gel(y,i)));
    2045      206784 :       gel(x,i) = rem(c, T);
    2046             :     }
    2047      288745 :     for (   ; i<=dx; i++)
    2048             :     {
    2049      198098 :       GEN c = gmul(y_lead, gel(x,i));
    2050      198098 :       gel(x,i) = rem(c, T);
    2051             :     }
    2052      100516 :     do { x++; dx--; } while (dx >= 0 && gequal0(gel(x,0)));
    2053       90647 :     if (dx < dy) break;
    2054       65703 :     if (gc_needed(av2,1))
    2055             :     {
    2056           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudorem dx = %ld >= %ld",dx,dy);
    2057           0 :       gerepilecoeffs(av2,x,dx+1);
    2058             :     }
    2059             :   }
    2060       24944 :   if (dx < 0) return pol_0(vx);
    2061       23411 :   lx = dx+3; x -= 2;
    2062       23411 :   x[0] = evaltyp(t_POL) | evallg(lx);
    2063       23411 :   x[1] = evalsigne(1) | evalvarn(vx);
    2064       23411 :   x = RgX_recip_shallow(x);
    2065       23411 :   if (p)
    2066             :   { /* multiply by y[0]^p   [beware dummy vars from FpX_FpXY_resultant] */
    2067        2730 :     GEN t = y_lead;
    2068        2730 :     if (T && typ(t) == t_POL && varn(t) == varn(T))
    2069           0 :       t = RgXQ_powu(t, p, T);
    2070             :     else
    2071        2730 :       t = gpowgs(t, p);
    2072        7496 :     for (i=2; i<lx; i++)
    2073             :     {
    2074        4766 :       GEN c = gmul(gel(x,i), t);
    2075        4766 :       gel(x,i) = rem(c,T);
    2076             :     }
    2077        2730 :     if (!T) return gerepileupto(av, x);
    2078             :   }
    2079       20681 :   return gerepilecopy(av, x);
    2080             : }
    2081             : 
    2082             : GEN
    2083       28122 : RgX_pseudorem(GEN x, GEN y) { return RgXQX_pseudorem(x,y, NULL); }
    2084             : 
    2085             : /* Compute z,r s.t lc(y)^(dx-dy+1) x = z y + r */
    2086             : GEN
    2087        7829 : RgXQX_pseudodivrem(GEN x, GEN y, GEN T, GEN *ptr)
    2088             : {
    2089        7829 :   long vx = varn(x), dx, dy, dz, i, iz, lx, lz, p;
    2090        7829 :   pari_sp av = avma, av2;
    2091             :   GEN z, r, ypow, y_lead;
    2092             : 
    2093        7829 :   if (!signe(y)) pari_err_INV("RgXQX_pseudodivrem",y);
    2094        7829 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2095        7829 :   if (gequal1(y_lead)) return T? RgXQX_divrem(x,y, T, ptr): RgX_divrem(x,y, ptr);
    2096        7257 :   dx = degpol(x);
    2097        7257 :   if (dx < dy) { *ptr = RgX_copy(x); return pol_0(vx); }
    2098        7257 :   if (dx == dy)
    2099             :   {
    2100          70 :     GEN x_lead = gel(x,lg(x)-1);
    2101          70 :     x = RgX_renormalize_lg(leafcopy(x), lg(x)-1);
    2102          70 :     y = RgX_renormalize_lg(leafcopy(y), lg(y)-1);
    2103          70 :     r = RgX_sub(RgX_Rg_mul(x, y_lead), RgX_Rg_mul(y, x_lead));
    2104          70 :     *ptr = gerepileupto(av, r); return scalarpol(x_lead, vx);
    2105             :   }
    2106        7187 :   (void)new_chunk(2);
    2107        7187 :   x = RgX_recip_shallow(x)+2;
    2108        7187 :   y = RgX_recip_shallow(y)+2;
    2109             :   /* pay attention to sparse divisors */
    2110       36794 :   for (i = 1; i <= dy; i++)
    2111       29607 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2112        7187 :   dz = dx-dy; p = dz+1;
    2113        7187 :   lz = dz+3;
    2114        7187 :   z = cgetg(lz, t_POL);
    2115        7187 :   z[1] = evalsigne(1) | evalvarn(vx);
    2116       27546 :   for (i = 2; i < lz; i++) gel(z,i) = gen_0;
    2117        7187 :   ypow = new_chunk(dz+1);
    2118        7187 :   gel(ypow,0) = gen_1;
    2119        7187 :   gel(ypow,1) = y_lead;
    2120       13172 :   for (i=2; i<=dz; i++)
    2121             :   {
    2122        5985 :     GEN c = gmul(gel(ypow,i-1), y_lead);
    2123        5985 :     gel(ypow,i) = rem(c,T);
    2124             :   }
    2125        7187 :   av2 = avma;
    2126        7187 :   for (iz=2;;)
    2127             :   {
    2128        7810 :     p--;
    2129       14997 :     gel(z,iz++) = rem(gmul(gel(x,0), gel(ypow,p)), T);
    2130       65083 :     for (i=1; i<=dy; i++)
    2131             :     {
    2132       50086 :       GEN c = gmul(y_lead, gel(x,i));
    2133       50086 :       if (gel(y,i)) c = gsub(c, gmul(gel(x,0),gel(y,i)));
    2134       50086 :       gel(x,i) = rem(c, T);
    2135             :     }
    2136       39138 :     for (   ; i<=dx; i++)
    2137             :     {
    2138       24141 :       GEN c = gmul(y_lead, gel(x,i));
    2139       24141 :       gel(x,i) = rem(c,T);
    2140             :     }
    2141       14997 :     x++; dx--;
    2142       20359 :     while (dx >= dy && gequal0(gel(x,0))) { x++; dx--; iz++; }
    2143       14997 :     if (dx < dy) break;
    2144        7810 :     if (gc_needed(av2,1))
    2145             :     {
    2146           0 :       GEN X = x-2;
    2147           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudodivrem dx=%ld >= %ld",dx,dy);
    2148           0 :       X[0] = evaltyp(t_POL)|evallg(dx+3); X[1] = z[1]; /* hack */
    2149           0 :       gerepileall(av2,2, &X, &z); x = X+2;
    2150             :     }
    2151             :   }
    2152       12941 :   while (dx >= 0 && gequal0(gel(x,0))) { x++; dx--; }
    2153        7187 :   if (dx < 0)
    2154          98 :     x = pol_0(vx);
    2155             :   else
    2156             :   {
    2157        7089 :     lx = dx+3; x -= 2;
    2158        7089 :     x[0] = evaltyp(t_POL) | evallg(lx);
    2159        7089 :     x[1] = evalsigne(1) | evalvarn(vx);
    2160        7089 :     x = RgX_recip_shallow(x);
    2161             :   }
    2162        7187 :   z = RgX_recip_shallow(z);
    2163        7187 :   r = x;
    2164        7187 :   if (p)
    2165             :   {
    2166        3234 :     GEN c = gel(ypow,p); r = RgX_Rg_mul(r, c);
    2167        3234 :     if (T && typ(c) == t_POL && varn(c) == varn(T)) r = RgXQX_red(r, T);
    2168             :   }
    2169        7187 :   gerepileall(av, 2, &z, &r);
    2170        7187 :   *ptr = r; return z;
    2171             : }
    2172             : GEN
    2173        7626 : RgX_pseudodivrem(GEN x, GEN y, GEN *ptr)
    2174        7626 : { return RgXQX_pseudodivrem(x,y,NULL,ptr); }
    2175             : 
    2176             : GEN
    2177           0 : RgXQX_mul(GEN x, GEN y, GEN T)
    2178             : {
    2179           0 :   return RgXQX_red(RgX_mul(x,y), T);
    2180             : }
    2181             : GEN
    2182    76878166 : RgX_Rg_mul(GEN y, GEN x) {
    2183             :   long i, ly;
    2184    76878166 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2185    76878160 :   if (ly == 2) return z;
    2186   294313326 :   for (i = 2; i < ly; i++) gel(z,i) = gmul(x,gel(y,i));
    2187    76791890 :   return normalizepol_lg(z,ly);
    2188             : }
    2189             : GEN
    2190         560 : RgX_muls(GEN y, long x) {
    2191             :   long i, ly;
    2192         560 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2193         560 :   if (ly == 2) return z;
    2194        3185 :   for (i = 2; i < ly; i++) gel(z,i) = gmulsg(x,gel(y,i));
    2195         462 :   return normalizepol_lg(z,ly);
    2196             : }
    2197             : GEN
    2198         756 : RgXQX_RgXQ_mul(GEN x, GEN y, GEN T)
    2199             : {
    2200         756 :   return RgXQX_red(RgX_Rg_mul(x,y), T);
    2201             : }
    2202             : GEN
    2203          77 : RgXQV_RgXQ_mul(GEN v, GEN x, GEN T)
    2204             : {
    2205          77 :   return RgXQV_red(RgV_Rg_mul(v,x), T);
    2206             : }
    2207             : 
    2208             : GEN
    2209           0 : RgXQX_sqr(GEN x, GEN T)
    2210             : {
    2211           0 :   return RgXQX_red(RgX_sqr(x), T);
    2212             : }
    2213             : 
    2214             : GEN
    2215           0 : RgXQX_powers(GEN P, long n, GEN T)
    2216             : {
    2217           0 :   GEN v = cgetg(n+2, t_VEC);
    2218             :   long i;
    2219           0 :   gel(v, 1) = pol_1(varn(T));
    2220           0 :   if (n==0) return v;
    2221           0 :   gel(v, 2) = gcopy(P);
    2222           0 :   for (i = 2; i <= n; i++) gel(v,i+1) = RgXQX_mul(P, gel(v,i), T);
    2223           0 :   return v;
    2224             : }
    2225             : 
    2226             : static GEN
    2227      142602 : _add(void *data, GEN x, GEN y) { (void)data; return RgX_add(x, y); }
    2228             : static GEN
    2229           0 : _sub(void *data, GEN x, GEN y) { (void)data; return RgX_sub(x, y); }
    2230             : static GEN
    2231       59561 : _sqr(void *data, GEN x) { return RgXQ_sqr(x, (GEN)data); }
    2232             : static GEN
    2233       49676 : _mul(void *data, GEN x, GEN y) { return RgXQ_mul(x,y, (GEN)data); }
    2234             : static GEN
    2235      255751 : _cmul(void *data, GEN P, long a, GEN x) { (void)data; return RgX_Rg_mul(x,gel(P,a+2)); }
    2236             : static GEN
    2237      255520 : _one(void *data) { return pol_1(varn((GEN)data)); }
    2238             : static GEN
    2239         469 : _zero(void *data) { return pol_0(varn((GEN)data)); }
    2240             : static GEN
    2241      154620 : _red(void *data, GEN x) { (void)data; return gcopy(x); }
    2242             : 
    2243             : static struct bb_algebra RgXQ_algebra = { _red, _add, _sub,
    2244             :               _mul, _sqr, _one, _zero };
    2245             : 
    2246             : GEN
    2247           0 : RgX_RgXQV_eval(GEN Q, GEN x, GEN T)
    2248             : {
    2249           0 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)T,&RgXQ_algebra,_cmul);
    2250             : }
    2251             : 
    2252             : GEN
    2253      112897 : RgX_RgXQ_eval(GEN Q, GEN x, GEN T)
    2254             : {
    2255      112897 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2256      112897 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)T,&RgXQ_algebra,_cmul);
    2257             : }
    2258             : 
    2259             : /* mod X^n */
    2260             : struct modXn {
    2261             :   long v; /* varn(X) */
    2262             :   long n;
    2263             : } ;
    2264             : static GEN
    2265       11613 : _sqrXn(void *data, GEN x) {
    2266       11613 :   struct modXn *S = (struct modXn*)data;
    2267       11613 :   return RgXn_sqr(x, S->n);
    2268             : }
    2269             : static GEN
    2270        4234 : _mulXn(void *data, GEN x, GEN y) {
    2271        4234 :   struct modXn *S = (struct modXn*)data;
    2272        4234 :   return RgXn_mul(x,y, S->n);
    2273             : }
    2274             : static GEN
    2275        1477 : _oneXn(void *data) {
    2276        1477 :   struct modXn *S = (struct modXn*)data;
    2277        1477 :   return pol_1(S->v);
    2278             : }
    2279             : static GEN
    2280           0 : _zeroXn(void *data) {
    2281           0 :   struct modXn *S = (struct modXn*)data;
    2282           0 :   return pol_0(S->v);
    2283             : }
    2284             : static struct bb_algebra RgXn_algebra = { _red, _add, _sub, _mulXn, _sqrXn,
    2285             :                                           _oneXn, _zeroXn };
    2286             : 
    2287             : GEN
    2288         357 : RgXn_powers(GEN x, long m, long n)
    2289             : {
    2290         357 :   long d = degpol(x);
    2291         357 :   int use_sqr = (d<<1) >= n;
    2292             :   struct modXn S;
    2293         357 :   S.v = varn(x); S.n = n;
    2294         357 :   return gen_powers(x,m,use_sqr,(void*)&S,_sqrXn,_mulXn,_oneXn);
    2295             : }
    2296             : 
    2297             : GEN
    2298        2216 : RgXn_powu_i(GEN x, ulong m, long n)
    2299             : {
    2300             :   struct modXn S;
    2301             :   long v;
    2302        2216 :   if (n == 1) return x;
    2303        2188 :   v = RgX_valrem(x, &x);
    2304        2188 :   if (v) n -= m * v;
    2305        2188 :   S.v = varn(x); S.n = n;
    2306        2188 :   x = gen_powu_i(x, m, (void*)&S,_sqrXn,_mulXn);
    2307        2188 :   if (v) x = RgX_shift_shallow(x, m * v);
    2308        2188 :   return x;
    2309             : }
    2310             : GEN
    2311           0 : RgXn_powu(GEN x, ulong m, long n)
    2312             : {
    2313             :   pari_sp av;
    2314           0 :   if (n == 1) return gcopy(x);
    2315           0 :   av = avma; return gerepilecopy(av, RgXn_powu_i(x, m, n));
    2316             : }
    2317             : 
    2318             : GEN
    2319         714 : RgX_RgXnV_eval(GEN Q, GEN x, long n)
    2320             : {
    2321             :   struct modXn S;
    2322         714 :   S.v = varn(gel(x,2)); S.n = n;
    2323         714 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)&S,&RgXn_algebra,_cmul);
    2324             : }
    2325             : 
    2326             : GEN
    2327           0 : RgX_RgXn_eval(GEN Q, GEN x, long n)
    2328             : {
    2329           0 :   int use_sqr = 2*degpol(x) >= n;
    2330             :   struct modXn S;
    2331           0 :   S.v = varn(x); S.n = n;
    2332           0 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2333             : }
    2334             : 
    2335             : /* Q(x) mod t^n, x in R[t], n >= 1 */
    2336             : GEN
    2337        3388 : RgXn_eval(GEN Q, GEN x, long n)
    2338             : {
    2339        3388 :   long d = degpol(x);
    2340             :   int use_sqr;
    2341             :   struct modXn S;
    2342        3388 :   if (d == 1 && isrationalzero(gel(x,2)))
    2343             :   {
    2344        3381 :     GEN y = RgX_unscale(Q, gel(x,3));
    2345        3381 :     setvarn(y, varn(x)); return y;
    2346             :   }
    2347           7 :   S.v = varn(x);
    2348           7 :   S.n = n;
    2349           7 :   use_sqr = (d<<1) >= n;
    2350           7 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2351             : }
    2352             : 
    2353             : /* (f*g mod t^n) \ t^n2, assuming 2*n2 >= n */
    2354             : static GEN
    2355     1475335 : RgXn_mulhigh(GEN f, GEN g, long n2, long n)
    2356             : {
    2357     1475335 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2358     1475335 :   return RgX_add(RgX_mulhigh_i(fl, g, n2), RgXn_mul(fh, g, n - n2));
    2359             : }
    2360             : 
    2361             : /* (f^2 mod t^n) \ t^n2, assuming 2*n2 >= n */
    2362             : static GEN
    2363           0 : RgXn_sqrhigh(GEN f, long n2, long n)
    2364             : {
    2365           0 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2366           0 :   return RgX_add(RgX_mulhigh_i(fl, f, n2), RgXn_mul(fh, f, n - n2));
    2367             : }
    2368             : 
    2369             : GEN
    2370      904975 : RgXn_inv_i(GEN f, long e)
    2371             : {
    2372             :   pari_sp av;
    2373             :   ulong mask;
    2374             :   GEN W, a;
    2375      904975 :   long v = varn(f), n = 1;
    2376             : 
    2377      904975 :   if (!signe(f)) pari_err_INV("RgXn_inv",f);
    2378      904968 :   a = ginv(gel(f,2));
    2379      904965 :   if (e == 1) return scalarpol(a, v);
    2380      901675 :   else if (e == 2)
    2381             :   {
    2382             :     GEN b;
    2383      214462 :     if (degpol(f) <= 0 || gequal0(b = gel(f,3))) return scalarpol(a, v);
    2384      192768 :     av = avma; b = gneg(b);
    2385      192769 :     if (!gequal1(a)) b = gmul(b, gsqr(a));
    2386      192767 :     W = deg1pol_shallow(b, a, v);
    2387      192767 :     return gcopy(W);
    2388             :   }
    2389      687213 :   W = scalarpol_shallow(a,v);
    2390      687213 :   mask = quadratic_prec_mask(e);
    2391      687213 :   av = avma;
    2392     2161309 :   for (;mask>1;)
    2393             :   {
    2394             :     GEN u, fr;
    2395     1474096 :     long n2 = n;
    2396     1474096 :     n<<=1; if (mask & 1) n--;
    2397     1474096 :     mask >>= 1;
    2398     1474096 :     fr = RgXn_red_shallow(f, n);
    2399     1474096 :     u = RgXn_mul(W, RgXn_mulhigh(fr, W, n2, n), n-n2);
    2400     1474096 :     W = RgX_sub(W, RgX_shift_shallow(u, n2));
    2401     1474096 :     if (gc_needed(av,2))
    2402             :     {
    2403           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_inv, e = %ld", n);
    2404           0 :       W = gerepileupto(av, W);
    2405             :     }
    2406             :   }
    2407      687213 :   return W;
    2408             : }
    2409             : 
    2410             : static GEN
    2411           0 : RgXn_inv_FpX(GEN x, long e, GEN p)
    2412             : {
    2413           0 :   pari_sp av = avma;
    2414             :   GEN r;
    2415           0 :   if (lgefint(p) == 3)
    2416             :   {
    2417           0 :     ulong pp = uel(p, 2);
    2418           0 :     if (pp == 2)
    2419           0 :       r = F2x_to_ZX(F2xn_inv(RgX_to_F2x(x), e));
    2420             :     else
    2421           0 :       r = Flx_to_ZX_inplace(Flxn_inv(RgX_to_Flx(x, pp), e, pp));
    2422             :   }
    2423             :   else
    2424           0 :     r = FpXn_inv(RgX_to_FpX(x, p), e, p);
    2425           0 :   return gerepileupto(av, FpX_to_mod(r, p));
    2426             : }
    2427             : 
    2428             : static GEN
    2429           0 : RgXn_inv_FpXQX(GEN x, long n, GEN pol, GEN p)
    2430             : {
    2431           0 :   pari_sp av = avma;
    2432           0 :   GEN r, T = RgX_to_FpX(pol, p);
    2433           0 :   if (signe(T) == 0) pari_err_OP("/", gen_1, x);
    2434           0 :   r = FpXQXn_inv(RgX_to_FpXQX(x, T, p), n, T, p);
    2435           0 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2436             : }
    2437             : 
    2438             : #define code(t1,t2) ((t1 << 6) | t2)
    2439             : 
    2440             : static GEN
    2441       77905 : RgXn_inv_fast(GEN x, long e)
    2442             : {
    2443             :   GEN p, pol;
    2444             :   long pa;
    2445       77905 :   long t = RgX_type(x,&p,&pol,&pa);
    2446       77905 :   switch(t)
    2447             :   {
    2448           0 :     case t_INTMOD: return RgXn_inv_FpX(x, e, p);
    2449           0 :     case code(t_POLMOD, t_INTMOD):
    2450           0 :                    return RgXn_inv_FpXQX(x, e, pol, p);
    2451       77905 :     default:       return NULL;
    2452             :   }
    2453             : }
    2454             : #undef code
    2455             : 
    2456             : GEN
    2457       77905 : RgXn_inv(GEN f, long e)
    2458             : {
    2459       77905 :   GEN h = RgXn_inv_fast(f, e);
    2460       77905 :   if (h) return h;
    2461       77905 :   return RgXn_inv_i(f, e);
    2462             : }
    2463             : 
    2464             : /* Compute intformal(x^n*S)/x^(n+1) */
    2465             : static GEN
    2466       14210 : RgX_integXn(GEN x, long n)
    2467             : {
    2468       14210 :   long i, lx = lg(x);
    2469             :   GEN y;
    2470       14210 :   if (lx == 2) return RgX_copy(x);
    2471       13846 :   y = cgetg(lx, t_POL); y[1] = x[1];
    2472       27713 :   for (i=2; i<lx; i++)
    2473       13867 :     gel(y,i) = gdivgs(gel(x,i), n+i-1);
    2474       13846 :   return RgX_renormalize_lg(y, lx);;
    2475             : }
    2476             : 
    2477             : GEN
    2478       12971 : RgXn_expint(GEN h, long e)
    2479             : {
    2480       12971 :   pari_sp av = avma, av2;
    2481       12971 :   long v = varn(h), n;
    2482       12971 :   GEN f = pol_1(v), g;
    2483             :   ulong mask;
    2484             : 
    2485       12971 :   if (!signe(h)) return f;
    2486       12971 :   g = pol_1(v);
    2487       12971 :   n = 1; mask = quadratic_prec_mask(e);
    2488       12971 :   av2 = avma;
    2489       14210 :   for (;mask>1;)
    2490             :   {
    2491             :     GEN u, w;
    2492       14210 :     long n2 = n;
    2493       14210 :     n<<=1; if (mask & 1) n--;
    2494       14210 :     mask >>= 1;
    2495       14210 :     u = RgXn_mul(g, RgX_mulhigh_i(f, RgXn_red_shallow(h, n2-1), n2-1), n-n2);
    2496       14210 :     u = RgX_add(u, RgX_shift_shallow(RgXn_red_shallow(h, n-1), 1-n2));
    2497       14210 :     w = RgXn_mul(f, RgX_integXn(u, n2-1), n-n2);
    2498       14210 :     f = RgX_add(f, RgX_shift_shallow(w, n2));
    2499       14210 :     if (mask<=1) break;
    2500        1239 :     u = RgXn_mul(g, RgXn_mulhigh(f, g, n2, n), n-n2);
    2501        1239 :     g = RgX_sub(g, RgX_shift_shallow(u, n2));
    2502        1239 :     if (gc_needed(av2,2))
    2503             :     {
    2504           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"RgXn_expint, e = %ld", n);
    2505           0 :       gerepileall(av2, 2, &f, &g);
    2506             :     }
    2507             :   }
    2508       12971 :   return gerepileupto(av, f);
    2509             : }
    2510             : 
    2511             : GEN
    2512           0 : RgXn_exp(GEN h, long e)
    2513             : {
    2514           0 :   long d = degpol(h);
    2515           0 :   if (d < 0) return pol_1(varn(h));
    2516           0 :   if (!d || !gequal0(gel(h,2)))
    2517           0 :     pari_err_DOMAIN("RgXn_exp","valuation", "<", gen_1, h);
    2518           0 :   return RgXn_expint(RgX_deriv(h), e);
    2519             : }
    2520             : 
    2521             : GEN
    2522          91 : RgXn_reverse(GEN f, long e)
    2523             : {
    2524          91 :   pari_sp av = avma, av2;
    2525             :   ulong mask;
    2526             :   GEN fi, a, df, W, an;
    2527          91 :   long v = varn(f), n=1;
    2528          91 :   if (degpol(f)<1 || !gequal0(gel(f,2)))
    2529           0 :     pari_err_INV("serreverse",f);
    2530          91 :   fi = ginv(gel(f,3));
    2531          91 :   a = deg1pol_shallow(fi,gen_0,v);
    2532          91 :   if (e <= 2) return gerepilecopy(av, a);
    2533          91 :   W = scalarpol(fi,v);
    2534          91 :   df = RgX_deriv(f);
    2535          91 :   mask = quadratic_prec_mask(e);
    2536          91 :   av2 = avma;
    2537         448 :   for (;mask>1;)
    2538             :   {
    2539             :     GEN u, fa, fr;
    2540         357 :     long n2 = n, rt;
    2541         357 :     n<<=1; if (mask & 1) n--;
    2542         357 :     mask >>= 1;
    2543         357 :     fr = RgXn_red_shallow(f, n);
    2544         357 :     rt = brent_kung_optpow(degpol(fr), 4, 3);
    2545         357 :     an = RgXn_powers(a, rt, n);
    2546         357 :     if (n>1)
    2547             :     {
    2548         357 :       long n4 = (n2+1)>>1;
    2549         357 :       GEN dfr = RgXn_red_shallow(df, n2);
    2550         357 :       dfr = RgX_RgXnV_eval(dfr, RgXnV_red_shallow(an, n2), n2);
    2551         357 :       u = RgX_shift(RgX_Rg_sub(RgXn_mul(W, dfr, n2), gen_1), -n4);
    2552         357 :       W = RgX_sub(W, RgX_shift(RgXn_mul(u, W, n2-n4), n4));
    2553             :     }
    2554         357 :     fa = RgX_sub(RgX_RgXnV_eval(fr, an, n), pol_x(v));
    2555         357 :     fa = RgX_shift(fa, -n2);
    2556         357 :     a = RgX_sub(a, RgX_shift(RgXn_mul(W, fa, n-n2), n2));
    2557         357 :     if (gc_needed(av2,2))
    2558             :     {
    2559           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_reverse, e = %ld", n);
    2560           0 :       gerepileall(av2, 2, &a, &W);
    2561             :     }
    2562             :   }
    2563          91 :   return gerepileupto(av, a);
    2564             : }
    2565             : 
    2566             : GEN
    2567           0 : RgXn_sqrt(GEN h, long e)
    2568             : {
    2569           0 :   pari_sp av = avma, av2;
    2570           0 :   long v = varn(h), n = 1;
    2571           0 :   GEN f = scalarpol(gen_1, v), df = f;
    2572           0 :   ulong mask = quadratic_prec_mask(e);
    2573           0 :   if (degpol(h)<0 || !gequal1(gel(h,2)))
    2574           0 :     pari_err_SQRTN("RgXn_sqrt",h);
    2575           0 :   av2 = avma;
    2576             :   while(1)
    2577           0 :   {
    2578           0 :     long n2 = n, m;
    2579             :     GEN g;
    2580           0 :     n<<=1; if (mask & 1) n--;
    2581           0 :     mask >>= 1;
    2582           0 :     m = n-n2;
    2583           0 :     g = RgX_sub(RgXn_sqrhigh(f, n2, n), RgX_shift_shallow(RgXn_red_shallow(h, n),-n2));
    2584           0 :     f = RgX_sub(f, RgX_shift_shallow(RgXn_mul(gmul2n(df, -1), g, m), n2));
    2585           0 :     if (mask==1) return gerepileupto(av, f);
    2586           0 :     g = RgXn_mul(df, RgXn_mulhigh(df, f, n2, n), m);
    2587           0 :     df = RgX_sub(df, RgX_shift_shallow(g, n2));
    2588           0 :     if (gc_needed(av2,2))
    2589             :     {
    2590           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_sqrt, e = %ld", n);
    2591           0 :       gerepileall(av2, 2, &f, &df);
    2592             :     }
    2593             :   }
    2594             : }
    2595             : 
    2596             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2597             : GEN
    2598       19140 : RgXQ_powu(GEN x, ulong n, GEN T)
    2599             : {
    2600       19140 :   pari_sp av = avma;
    2601             : 
    2602       19140 :   if (!n) return pol_1(varn(x));
    2603       17277 :   if (n == 1) return RgX_copy(x);
    2604       15302 :   x = gen_powu_i(x, n, (void*)T, &_sqr, &_mul);
    2605       15302 :   return gerepilecopy(av, x);
    2606             : }
    2607             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2608             : GEN
    2609       20837 : RgXQ_pow(GEN x, GEN n, GEN T)
    2610             : {
    2611             :   pari_sp av;
    2612       20837 :   long s = signe(n);
    2613             : 
    2614       20837 :   if (!s) return pol_1(varn(x));
    2615       20837 :   if (is_pm1(n) == 1)
    2616           0 :     return (s < 0)? RgXQ_inv(x, T): RgX_copy(x);
    2617       20837 :   av = avma;
    2618       20837 :   if (s < 0) x = RgXQ_inv(x, T);
    2619       20837 :   x = gen_pow_i(x, n, (void*)T, &_sqr, &_mul);
    2620       20837 :   return gerepilecopy(av, x);
    2621             : }
    2622             : static GEN
    2623      161481 : _ZXQsqr(void *data, GEN x) { return ZXQ_sqr(x, (GEN)data); }
    2624             : static GEN
    2625       79959 : _ZXQmul(void *data, GEN x, GEN y) { return ZXQ_mul(x,y, (GEN)data); }
    2626             : 
    2627             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2628             : GEN
    2629        2429 : ZXQ_powers(GEN x, long l, GEN T)
    2630             : {
    2631        2429 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2632        2429 :   return gen_powers(x, l, use_sqr, (void *)T,_ZXQsqr,_ZXQmul,_one);
    2633             : }
    2634             : 
    2635             : /* x,T in Z[X], n in N, compute lift(x^n mod T)) */
    2636             : GEN
    2637      163655 : ZXQ_powu(GEN x, ulong n, GEN T)
    2638             : {
    2639      163655 :   pari_sp av = avma;
    2640             : 
    2641      163655 :   if (!n) return pol_1(varn(x));
    2642      163655 :   if (n == 1) return ZX_copy(x);
    2643      107403 :   x = gen_powu_i(x, n, (void*)T, &_ZXQsqr, &_ZXQmul);
    2644      107411 :   return gerepilecopy(av, x);
    2645             : }
    2646             : 
    2647             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2648             : GEN
    2649        5369 : RgXQ_powers(GEN x, long l, GEN T)
    2650             : {
    2651        5369 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2652        5369 :   return gen_powers(x, l, use_sqr, (void *)T,_sqr,_mul,_one);
    2653             : }
    2654             : 
    2655             : /* a in K = Q[X]/(T), returns [a^0, ..., a^n] */
    2656             : GEN
    2657        3500 : QXQ_powers(GEN a, long n, GEN T)
    2658             : {
    2659        3500 :   GEN den, v = RgXQ_powers(Q_remove_denom(a, &den), n, T);
    2660             :   /* den*a integral; v[i+1] = (den*a)^i in K */
    2661        3500 :   if (den)
    2662             :   { /* restore denominators */
    2663        2156 :     GEN d = den;
    2664             :     long i;
    2665        2156 :     gel(v,2) = a;
    2666        6909 :     for (i=3; i<=n+1; i++) {
    2667        4753 :       d = mulii(d,den);
    2668        4753 :       gel(v,i) = RgX_Rg_div(gel(v,i), d);
    2669             :     }
    2670             :   }
    2671        3500 :   return v;
    2672             : }
    2673             : 
    2674             : static GEN
    2675        2016 : do_QXQ_eval(GEN v, long imin, GEN a, GEN T)
    2676             : {
    2677        2016 :   long l, i, m = 0;
    2678             :   GEN dz, z;
    2679        2016 :   GEN V = cgetg_copy(v, &l);
    2680        6496 :   for (i = imin; i < l; i++)
    2681             :   {
    2682        4480 :     GEN c = gel(v, i);
    2683        4480 :     if (typ(c) == t_POL) m = maxss(m, degpol(c));
    2684             :   }
    2685        2016 :   z = Q_remove_denom(QXQ_powers(a, m, T), &dz);
    2686        2079 :   for (i = 1; i < imin; i++) V[i] = v[i];
    2687        6496 :   for (i = imin; i < l; i++)
    2688             :   {
    2689        4480 :     GEN c = gel(v,i);
    2690        4480 :     if (typ(c) == t_POL) c = QX_ZXQV_eval(c, z, dz);
    2691        4480 :     gel(V,i) = c;
    2692             :   }
    2693        2016 :   return V;
    2694             : }
    2695             : /* [ s(a mod T) | s <- lift(v) ], a,T are QX, v a QXV */
    2696             : GEN
    2697        1953 : QXV_QXQ_eval(GEN v, GEN a, GEN T)
    2698        1953 : { return do_QXQ_eval(v, 1, a, T); }
    2699             : GEN
    2700          63 : QXX_QXQ_eval(GEN v, GEN a, GEN T)
    2701          63 : { return normalizepol(do_QXQ_eval(v, 2, a, T)); }
    2702             : 
    2703             : GEN
    2704         609 : RgXQ_matrix_pow(GEN y, long n, long m, GEN P)
    2705             : {
    2706         609 :   return RgXV_to_RgM(RgXQ_powers(y,m-1,P),n);
    2707             : }
    2708             : 
    2709             : GEN
    2710        5028 : RgXQ_norm(GEN x, GEN T)
    2711             : {
    2712             :   pari_sp av;
    2713        5028 :   long dx = degpol(x);
    2714             :   GEN L, y;
    2715             : 
    2716        5028 :   av = avma; y = resultant(T, x);
    2717        5028 :   L = leading_coeff(T);
    2718        5028 :   if (gequal1(L) || !signe(x)) return y;
    2719           0 :   return gerepileupto(av, gdiv(y, gpowgs(L, dx)));
    2720             : }
    2721             : 
    2722             : GEN
    2723     1622385 : RgX_blocks(GEN P, long n, long m)
    2724             : {
    2725     1622385 :   GEN z = cgetg(m+1,t_VEC);
    2726     1622385 :   long i,j, k=2, l = lg(P);
    2727     5060327 :   for(i=1; i<=m; i++)
    2728             :   {
    2729     3437942 :     GEN zi = cgetg(n+2,t_POL);
    2730     3437942 :     zi[1] = P[1];
    2731     3437942 :     gel(z,i) = zi;
    2732    11164245 :     for(j=2; j<n+2; j++)
    2733     7726303 :       gel(zi, j) = k==l ? gen_0 : gel(P,k++);
    2734     3437942 :     zi = RgX_renormalize_lg(zi, n+2);
    2735             :   }
    2736     1622385 :   return z;
    2737             : }
    2738             : 
    2739             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
    2740             : void
    2741     6212839 : RgX_even_odd(GEN p, GEN *pe, GEN *po)
    2742             : {
    2743     6212839 :   long n = degpol(p), v = varn(p), n0, n1, i;
    2744             :   GEN p0, p1;
    2745             : 
    2746     6212844 :   if (n <= 0) { *pe = RgX_copy(p); *po = zeropol(v); return; }
    2747             : 
    2748     6212557 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
    2749     6212557 :   p0 = cgetg(n0+2, t_POL); p0[1] = evalvarn(v)|evalsigne(1);
    2750     6212557 :   p1 = cgetg(n1+2, t_POL); p1[1] = evalvarn(v)|evalsigne(1);
    2751    25772972 :   for (i=0; i<n1; i++)
    2752             :   {
    2753    19560414 :     p0[2+i] = p[2+(i<<1)];
    2754    19560414 :     p1[2+i] = p[3+(i<<1)];
    2755             :   }
    2756     6212558 :   if (n1 != n0)
    2757     3058540 :     p0[2+i] = p[2+(i<<1)];
    2758     6212558 :   *pe = normalizepol(p0);
    2759     6212559 :   *po = normalizepol(p1);
    2760             : }
    2761             : 
    2762             : /* write p(X) = a_0(X^k) + Xa_1(X^k) + ... + X^(k-1)a_{k-1}(X^k), shallow function */
    2763             : GEN
    2764       40922 : RgX_splitting(GEN p, long k)
    2765             : {
    2766       40922 :   long n = degpol(p), v = varn(p), m, i, j, l;
    2767             :   GEN r;
    2768             : 
    2769       40922 :   m = n/k;
    2770       40922 :   r = cgetg(k+1,t_VEC);
    2771      225666 :   for(i=1; i<=k; i++)
    2772             :   {
    2773      184744 :     gel(r,i) = cgetg(m+3, t_POL);
    2774      184744 :     mael(r,i,1) = evalvarn(v)|evalsigne(1);
    2775             :   }
    2776      554981 :   for (j=1, i=0, l=2; i<=n; i++)
    2777             :   {
    2778      514059 :     gmael(r,j,l) = gel(p,2+i);
    2779      514059 :     if (j==k) { j=1; l++; } else j++;
    2780             :   }
    2781      225666 :   for(i=1; i<=k; i++)
    2782      184744 :     gel(r,i) = normalizepol_lg(gel(r,i),i<j?l+1:l);
    2783       40922 :   return r;
    2784             : }
    2785             : 
    2786             : /*******************************************************************/
    2787             : /*                                                                 */
    2788             : /*                        Kronecker form                           */
    2789             : /*                                                                 */
    2790             : /*******************************************************************/
    2791             : 
    2792             : /* z in R[Y] representing an elt in R[X,Y] mod T(Y) in Kronecker form,
    2793             :  * i.e subst(lift(z), x, y^(2deg(z)-1)). Recover the "real" z, with
    2794             :  * normalized coefficients */
    2795             : GEN
    2796      174925 : Kronecker_to_mod(GEN z, GEN T)
    2797             : {
    2798      174925 :   long i,j,lx,l = lg(z), N = (degpol(T)<<1) + 1;
    2799      174925 :   GEN x, t = cgetg(N,t_POL);
    2800      174925 :   t[1] = T[1];
    2801      174925 :   lx = (l-2) / (N-2); x = cgetg(lx+3,t_POL);
    2802      174925 :   x[1] = z[1];
    2803      174925 :   T = RgX_copy(T);
    2804     1936201 :   for (i=2; i<lx+2; i++, z+= N-2)
    2805             :   {
    2806     7938192 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
    2807     1761276 :     gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2808             :   }
    2809      174925 :   N = (l-2) % (N-2) + 2;
    2810      453055 :   for (j=2; j<N; j++) t[j] = z[j];
    2811      174925 :   gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2812      174925 :   return normalizepol_lg(x, i+1);
    2813             : }
    2814             : 
    2815             : /*******************************************************************/
    2816             : /*                                                                 */
    2817             : /*                        Domain detection                         */
    2818             : /*                                                                 */
    2819             : /*******************************************************************/
    2820             : 
    2821             : static GEN
    2822      151005 : zero_FpX_mod(GEN p, long v)
    2823             : {
    2824      151005 :   GEN r = cgetg(3,t_POL);
    2825      151005 :   r[1] = evalvarn(v);
    2826      151005 :   gel(r,2) = mkintmod(gen_0, icopy(p));
    2827      151005 :   return r;
    2828             : }
    2829             : 
    2830             : static GEN
    2831      499694 : RgX_mul_FpX(GEN x, GEN y, GEN p)
    2832             : {
    2833      499694 :   pari_sp av = avma;
    2834             :   GEN r;
    2835      499694 :   if (lgefint(p) == 3)
    2836             :   {
    2837      469562 :     ulong pp = uel(p, 2);
    2838      469562 :     r = Flx_to_ZX_inplace(Flx_mul(RgX_to_Flx(x, pp),
    2839             :                                   RgX_to_Flx(y, pp), pp));
    2840             :   }
    2841             :   else
    2842       30132 :     r = FpX_mul(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2843      499694 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2844      386377 :   return gerepileupto(av, FpX_to_mod(r, p));
    2845             : }
    2846             : 
    2847             : static GEN
    2848         532 : zero_FpXQX_mod(GEN pol, GEN p, long v)
    2849             : {
    2850         532 :   GEN r = cgetg(3,t_POL);
    2851         532 :   r[1] = evalvarn(v);
    2852         532 :   gel(r,2) = mkpolmod(mkintmod(gen_0, icopy(p)), gcopy(pol));
    2853         532 :   return r;
    2854             : }
    2855             : 
    2856             : static GEN
    2857        2037 : RgX_mul_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2858             : {
    2859        2037 :   pari_sp av = avma;
    2860             :   long dT;
    2861             :   GEN kx, ky, r;
    2862        2037 :   GEN T = RgX_to_FpX(pol, p);
    2863        2037 :   if (signe(T)==0) pari_err_OP("*", x, y);
    2864        2030 :   dT = degpol(T);
    2865        2030 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2866        2030 :   ky = ZXX_to_Kronecker(RgX_to_FpXQX(y, T, p), dT);
    2867        2030 :   r = FpX_mul(kx, ky, p);
    2868        2030 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2869        1505 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2870             : }
    2871             : 
    2872             : static GEN
    2873      404988 : RgX_liftred(GEN x, GEN T)
    2874      404988 : { return RgXQX_red(liftpol_shallow(x), T); }
    2875             : 
    2876             : static GEN
    2877      159853 : RgX_mul_QXQX(GEN x, GEN y, GEN T)
    2878             : {
    2879      159853 :   pari_sp av = avma;
    2880      159853 :   long dT = degpol(T);
    2881      159853 :   GEN r = QX_mul(ZXX_to_Kronecker(RgX_liftred(x, T), dT),
    2882             :                  ZXX_to_Kronecker(RgX_liftred(y, T), dT));
    2883      159853 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2884             : }
    2885             : 
    2886             : static GEN
    2887        1330 : RgX_sqr_FpX(GEN x, GEN p)
    2888             : {
    2889        1330 :   pari_sp av = avma;
    2890             :   GEN r;
    2891        1330 :   if (lgefint(p) == 3)
    2892             :   {
    2893        1126 :     ulong pp = uel(p, 2);
    2894        1126 :     r = Flx_to_ZX_inplace(Flx_sqr(RgX_to_Flx(x, pp), pp));
    2895             :   }
    2896             :   else
    2897         204 :     r = FpX_sqr(RgX_to_FpX(x, p), p);
    2898        1330 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2899        1232 :   return gerepileupto(av, FpX_to_mod(r, p));
    2900             : }
    2901             : 
    2902             : static GEN
    2903         196 : RgX_sqr_FpXQX(GEN x, GEN pol, GEN p)
    2904             : {
    2905         196 :   pari_sp av = avma;
    2906             :   long dT;
    2907         196 :   GEN kx, r, T = RgX_to_FpX(pol, p);
    2908         196 :   if (signe(T)==0) pari_err_OP("*",x,x);
    2909         189 :   dT = degpol(T);
    2910         189 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2911         189 :   r = FpX_sqr(kx, p);
    2912         189 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2913         189 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2914             : }
    2915             : 
    2916             : static GEN
    2917       13378 : RgX_sqr_QXQX(GEN x, GEN T)
    2918             : {
    2919       13378 :   pari_sp av = avma;
    2920       13378 :   long dT = degpol(T);
    2921       13378 :   GEN r = QX_sqr(ZXX_to_Kronecker(RgX_liftred(x, T), dT));
    2922       13378 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2923             : }
    2924             : 
    2925             : static GEN
    2926       67466 : RgX_rem_FpX(GEN x, GEN y, GEN p)
    2927             : {
    2928       67466 :   pari_sp av = avma;
    2929             :   GEN r;
    2930       67466 :   if (lgefint(p) == 3)
    2931             :   {
    2932       50676 :     ulong pp = uel(p, 2);
    2933       50676 :     r = Flx_to_ZX_inplace(Flx_rem(RgX_to_Flx(x, pp),
    2934             :                                   RgX_to_Flx(y, pp), pp));
    2935             :   }
    2936             :   else
    2937       16790 :     r = FpX_rem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2938       67466 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2939       29876 :   return gerepileupto(av, FpX_to_mod(r, p));
    2940             : }
    2941             : 
    2942             : static GEN
    2943       35952 : RgX_rem_QXQX(GEN x, GEN y, GEN T)
    2944             : {
    2945       35952 :   pari_sp av = avma;
    2946             :   GEN r;
    2947       35952 :   r = RgXQX_rem(RgX_liftred(x, T), RgX_liftred(y, T), T);
    2948       35952 :   return gerepilecopy(av, QXQX_to_mod_shallow(r, T));
    2949             : }
    2950             : static GEN
    2951          70 : RgX_rem_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2952             : {
    2953          70 :   pari_sp av = avma;
    2954             :   GEN r;
    2955          70 :   GEN T = RgX_to_FpX(pol, p);
    2956          70 :   if (signe(T) == 0) pari_err_OP("%", x, y);
    2957          63 :   if (lgefint(p) == 3)
    2958             :   {
    2959          55 :     ulong pp = uel(p, 2);
    2960          55 :     GEN Tp = ZX_to_Flx(T, pp);
    2961          55 :     r = FlxX_to_ZXX(FlxqX_rem(RgX_to_FlxqX(x, Tp, pp),
    2962             :                               RgX_to_FlxqX(y, Tp, pp), Tp, pp));
    2963             :   }
    2964             :   else
    2965           8 :     r = FpXQX_rem(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), T, p);
    2966          63 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2967          56 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2968             : }
    2969             : 
    2970             : #define code(t1,t2) ((t1 << 6) | t2)
    2971             : static GEN
    2972    53288257 : RgX_mul_fast(GEN x, GEN y)
    2973             : {
    2974             :   GEN p, pol;
    2975             :   long pa;
    2976    53288257 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2977    53288257 :   switch(t)
    2978             :   {
    2979    34504112 :     case t_INT:    return ZX_mul(x,y);
    2980     2384485 :     case t_FRAC:   return QX_mul(x,y);
    2981      104693 :     case t_FFELT:  return FFX_mul(x, y, pol);
    2982      499694 :     case t_INTMOD: return RgX_mul_FpX(x, y, p);
    2983      159853 :     case code(t_POLMOD, t_INT):
    2984             :     case code(t_POLMOD, t_FRAC):
    2985      159853 :                    return RgX_mul_QXQX(x, y, pol);
    2986        2037 :     case code(t_POLMOD, t_INTMOD):
    2987        2037 :                    return RgX_mul_FpXQX(x, y, pol, p);
    2988    15633383 :     default:       return NULL;
    2989             :   }
    2990             : }
    2991             : static GEN
    2992      142160 : RgX_sqr_fast(GEN x)
    2993             : {
    2994             :   GEN p, pol;
    2995             :   long pa;
    2996      142160 :   long t = RgX_type(x,&p,&pol,&pa);
    2997      142160 :   switch(t)
    2998             :   {
    2999      103883 :     case t_INT:    return ZX_sqr(x);
    3000       11599 :     case t_FRAC:   return QX_sqr(x);
    3001        2296 :     case t_FFELT:  return FFX_sqr(x, pol);
    3002        1330 :     case t_INTMOD: return RgX_sqr_FpX(x, p);
    3003       13378 :     case code(t_POLMOD, t_INT):
    3004             :     case code(t_POLMOD, t_FRAC):
    3005       13378 :                    return RgX_sqr_QXQX(x, pol);
    3006         196 :     case code(t_POLMOD, t_INTMOD):
    3007         196 :                    return RgX_sqr_FpXQX(x, pol, p);
    3008        9478 :     default:       return NULL;
    3009             :   }
    3010             : }
    3011             : 
    3012             : static GEN
    3013     7465642 : RgX_rem_fast(GEN x, GEN y)
    3014             : {
    3015             :   GEN p, pol;
    3016             :   long pa;
    3017     7465642 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    3018     7465642 :   switch(t)
    3019             :   {
    3020     4432965 :     case t_INT:    return ZX_is_monic(y) ? ZX_rem(x,y): NULL;
    3021     1542396 :     case t_FRAC:   return RgX_is_ZX(y) && ZX_is_monic(y) ? QX_ZX_rem(x,y): NULL;
    3022          84 :     case t_FFELT:  return FFX_rem(x, y, pol);
    3023       67466 :     case t_INTMOD: return RgX_rem_FpX(x, y, p);
    3024       35952 :     case code(t_POLMOD, t_INT):
    3025             :     case code(t_POLMOD, t_FRAC):
    3026       35952 :                    return RgX_rem_QXQX(x, y, pol);
    3027          70 :     case code(t_POLMOD, t_INTMOD):
    3028          70 :                    return RgX_rem_FpXQX(x, y, pol, p);
    3029     1386709 :     default:       return NULL;
    3030             :   }
    3031             : }
    3032             : 
    3033             : #undef code
    3034             : 
    3035             : GEN
    3036    45394301 : RgX_mul(GEN x, GEN y)
    3037             : {
    3038    45394301 :   GEN z = RgX_mul_fast(x,y);
    3039    45394289 :   if (!z) z = RgX_mul_i(x,y);
    3040    45394289 :   return z;
    3041             : }
    3042             : 
    3043             : GEN
    3044      130302 : RgX_sqr(GEN x)
    3045             : {
    3046      130302 :   GEN z = RgX_sqr_fast(x);
    3047      130295 :   if (!z) z = RgX_sqr_i(x);
    3048      130295 :   return z;
    3049             : }
    3050             : 
    3051             : GEN
    3052     7465642 : RgX_rem(GEN x, GEN y)
    3053             : {
    3054     7465642 :   GEN z = RgX_rem_fast(x, y);
    3055     7465635 :   if (!z) z = RgX_divrem_i(x, y, ONLY_REM);
    3056     7465635 :   return z;
    3057             : }
    3058             : 
    3059             : GEN
    3060     6404411 : RgXn_mul(GEN f, GEN g, long n)
    3061             : {
    3062     6404411 :   pari_sp av = avma;
    3063     6404411 :   GEN h = RgX_mul_fast(f,g);
    3064     6404411 :   if (!h) return RgXn_mul2(f,g,n);
    3065     1252138 :   if (degpol(h) < n) return h;
    3066      480864 :   return gerepilecopy(av, RgXn_red_shallow(h, n));
    3067             : }
    3068             : 
    3069             : GEN
    3070       11858 : RgXn_sqr(GEN f, long n)
    3071             : {
    3072       11858 :   pari_sp av = avma;
    3073       11858 :   GEN g = RgX_sqr_fast(f);
    3074       11858 :   if (!g) return RgXn_sqr2(f,n);
    3075       11613 :   if (degpol(g) < n) return g;
    3076       10521 :   return gerepilecopy(av, RgXn_red_shallow(g, n));
    3077             : }
    3078             : 
    3079             : /* (f*g) \/ x^n */
    3080             : GEN
    3081     1489545 : RgX_mulhigh_i(GEN f, GEN g, long n)
    3082             : {
    3083     1489545 :   GEN h = RgX_mul_fast(f,g);
    3084     1489545 :   return h? RgX_shift_shallow(h, -n): RgX_mulhigh_i2(f,g,n);
    3085             : }
    3086             : 
    3087             : /* (f*g) \/ x^n */
    3088             : GEN
    3089           0 : RgX_sqrhigh_i(GEN f, long n)
    3090             : {
    3091           0 :   GEN h = RgX_sqr_fast(f);
    3092           0 :   return h? RgX_shift_shallow(h, -n): RgX_sqrhigh_i2(f,n);
    3093             : }

Generated by: LCOV version 1.13