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 24038-ebe36f6c4) Lines: 1477 1643 89.9 %
Date: 2019-07-23 05:53:17 Functions: 178 193 92.2 %
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     3808876 : brent_kung_optpow(long d, long n, long m)
      29             : {
      30             :   long p, r;
      31     3808876 :   long pold=1, rold=n*(d-1);
      32    24276349 :   for(p=2; p<=d; p++)
      33             :   {
      34    20467473 :     r = m*(p-1) + n*((d-1)/p);
      35    20467473 :     if (r<rold) { pold=p; rold=r; }
      36             :   }
      37     3808876 :   return pold;
      38             : }
      39             : 
      40             : static GEN
      41     7035688 : 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     7035688 :   pari_sp av = avma;
      45             :   long i;
      46     7035688 :   GEN z = cmul(E,P,a,ff->one(E));
      47     7036536 :   if (!z) z = gen_0;
      48    49575092 :   for (i=1; i<=n; i++)
      49             :   {
      50    42539454 :     GEN t = cmul(E,P,a+i,gel(V,i+1));
      51    42544865 :     if (t) {
      52    31600732 :       z = ff->add(E, z, t);
      53    31594341 :       if (gc_needed(av,2)) z = gerepileupto(av, z);
      54             :     }
      55             :   }
      56     7035638 :   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     5891661 : 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     5891661 :   pari_sp av = avma;
      69     5891661 :   long l = lg(V)-1;
      70             :   GEN z, u;
      71             : 
      72     5891661 :   if (d < 0) return ff->zero(E);
      73     5573058 :   if (d < l) return gerepileupto(av, gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul));
      74      651665 :   if (l<2) pari_err_DOMAIN("gen_RgX_bkeval_powers", "#powers", "<",gen_2,V);
      75      651665 :   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      651665 :   d -= l;
      81      651665 :   z = gen_RgXQ_eval_powers(P,V,d+1,l-1,E,ff,cmul);
      82     2114364 :   while (d >= l-1)
      83             :   {
      84      811042 :     d -= l-1;
      85      811042 :     u = gen_RgXQ_eval_powers(P,V,d+1,l-2,E,ff,cmul);
      86      810991 :     z = ff->add(E,u, ff->mul(E,z,gel(V,l)));
      87      811016 :     if (gc_needed(av,2))
      88          92 :       z = gerepileupto(av, z);
      89             :   }
      90      651664 :   u = gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul);
      91      651659 :   z = ff->add(E,u, ff->mul(E,z,gel(V,d+2)));
      92      651662 :   return gerepileupto(av, ff->red(E,z));
      93             : }
      94             : 
      95             : GEN
      96      405734 : 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      405734 :   pari_sp av = avma;
     100             :   GEN z, V;
     101             :   long rtd;
     102      405734 :   if (d < 0) return ff->zero(E);
     103      405076 :   rtd = (long) sqrt((double)d);
     104      405076 :   V = gen_powers(x,rtd,use_sqr,E,ff->sqr,ff->mul,ff->one);
     105      405074 :   z = gen_bkeval_powers(Q, d, V, E, ff, cmul);
     106      405075 :   return gerepileupto(av, z);
     107             : }
     108             : 
     109             : static GEN
     110      928812 : _gen_nored(void *E, GEN x) { (void)E; return x; }
     111             : static GEN
     112    12841711 : _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      782140 : _gen_mul(void *E, GEN x, GEN y) { (void)E; return gmul(x, y); }
     117             : static GEN
     118      256622 : _gen_sqr(void *E, GEN x) { (void)E; return gsqr(x); }
     119             : static GEN
     120      947194 : _gen_one(void *E) { (void)E; return gen_1; }
     121             : static GEN
     122       21035 : _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      476553 : _gen_cmul(void *E, GEN P, long a, GEN x)
     129      476553 : {(void)E; return gmul(gel(P,a+2), x);}
     130             : 
     131             : GEN
     132      154406 : RgX_RgV_eval(GEN Q, GEN x)
     133             : {
     134      154406 :   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        2905 : RgXV_RgV_eval(GEN Q, GEN x)
     145             : {
     146        2905 :   long i, l = lg(Q), vQ = gvar(Q);
     147        2905 :   GEN v = cgetg(l, t_VEC);
     148      243236 :   for (i = 1; i < l; i++)
     149             :   {
     150      240331 :     GEN Qi = gel(Q, i);
     151      240331 :     gel(v, i) = typ(Qi)==t_POL && varn(Qi)==vQ? RgX_RgV_eval(Qi, x): gcopy(Qi);
     152             :   }
     153        2905 :   return v;
     154             : }
     155             : 
     156             : const struct bb_algebra *
     157      118402 : get_Rg_algebra(void)
     158             : {
     159      118402 :   return &Rg_algebra;
     160             : }
     161             : 
     162             : static struct bb_ring Rg_ring = {  _gen_add, _gen_mul, _gen_sqr };
     163             : 
     164             : static GEN
     165       17766 : _RgX_divrem(void *E, GEN x, GEN y, GEN *r)
     166             : {
     167             :   (void) E;
     168       17766 :   return RgX_divrem(x, y, r);
     169             : }
     170             : 
     171             : GEN
     172        4354 : RgX_digits(GEN x, GEN T)
     173             : {
     174        4354 :   pari_sp av = avma;
     175        4354 :   long d = degpol(T), n = (lgpol(x)+d-1)/d;
     176        4354 :   GEN z = gen_digits(x,T,n,NULL, &Rg_ring, _RgX_divrem);
     177        4354 :   return gerepileupto(av, z);
     178             : }
     179             : 
     180             : /*******************************************************************/
     181             : /*                                                                 */
     182             : /*                         RgX                                     */
     183             : /*                                                                 */
     184             : /*******************************************************************/
     185             : 
     186             : long
     187    24567034 : RgX_equal(GEN x, GEN y)
     188             : {
     189    24567034 :   long i = lg(x);
     190             : 
     191    24567034 :   if (i != lg(y)) return 0;
     192   103523227 :   for (i--; i > 1; i--)
     193    79084623 :     if (!gequal(gel(x,i),gel(y,i))) return 0;
     194    24438604 :   return 1;
     195             : }
     196             : 
     197             : /* Returns 1 in the base ring over which x is defined */
     198             : /* HACK: this also works for t_SER */
     199             : GEN
     200     2889285 : Rg_get_1(GEN x)
     201             : {
     202             :   GEN p, T;
     203     2889285 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     204     2889285 :   if (RgX_type_is_composite(tx))
     205       13545 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     206     2889285 :   switch(tx)
     207             :   {
     208         301 :     case t_INTMOD: retmkintmod(is_pm1(p)? gen_0: gen_1, icopy(p));
     209           7 :     case t_PADIC: return cvtop(gen_1, p, lx);
     210         735 :     case t_FFELT: return FF_1(T);
     211     2888242 :     default: return gen_1;
     212             :   }
     213             : }
     214             : /* Returns 0 in the base ring over which x is defined */
     215             : /* HACK: this also works for t_SER */
     216             : GEN
     217     6971741 : Rg_get_0(GEN x)
     218             : {
     219             :   GEN p, T;
     220     6971741 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     221     6971741 :   if (RgX_type_is_composite(tx))
     222       50554 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     223     6971741 :   switch(tx)
     224             :   {
     225         322 :     case t_INTMOD: retmkintmod(gen_0, icopy(p));
     226           0 :     case t_PADIC: return cvtop(gen_0, p, lx);
     227         231 :     case t_FFELT: return FF_zero(T);
     228     6971188 :     default: return gen_0;
     229             :   }
     230             : }
     231             : 
     232             : GEN
     233        4291 : QX_ZXQV_eval(GEN P, GEN V, GEN dV)
     234             : {
     235        4291 :   long i, n = degpol(P);
     236             :   GEN z, dz, dP;
     237        4291 :   if (n < 0) return gen_0;
     238        4291 :   P = Q_remove_denom(P, &dP);
     239        4291 :   z = gel(P,2); if (n == 0) return icopy(z);
     240        2345 :   if (dV) z = mulii(dV, z); /* V[1] = dV */
     241        2345 :   z = ZX_Z_add_shallow(ZX_Z_mul(gel(V,2),gel(P,3)), z);
     242        2345 :   for (i=2; i<=n; i++) z = ZX_add(ZX_Z_mul(gel(V,i+1),gel(P,2+i)), z);
     243        2345 :   dz = mul_denom(dP, dV);
     244        2345 :   return dz? RgX_Rg_div(z, dz): z;
     245             : }
     246             : 
     247             : /* Return P(h * x), not memory clean */
     248             : GEN
     249       11067 : RgX_unscale(GEN P, GEN h)
     250             : {
     251       11067 :   long i, l = lg(P);
     252       11067 :   GEN hi = gen_1, Q = cgetg(l, t_POL);
     253       11067 :   Q[1] = P[1];
     254       11067 :   if (l == 2) return Q;
     255       11039 :   gel(Q,2) = gcopy(gel(P,2));
     256       44982 :   for (i=3; i<l; i++)
     257             :   {
     258       33943 :     hi = gmul(hi,h);
     259       33943 :     gel(Q,i) = gmul(gel(P,i), hi);
     260             :   }
     261       11039 :   return Q;
     262             : }
     263             : /* P a ZX, Return P(h * x), not memory clean; optimize for h = -1 */
     264             : GEN
     265      371902 : ZX_z_unscale(GEN P, long h)
     266             : {
     267      371902 :   long i, l = lg(P);
     268      371902 :   GEN Q = cgetg(l, t_POL);
     269      371902 :   Q[1] = P[1];
     270      371902 :   if (l == 2) return Q;
     271      371902 :   gel(Q,2) = gel(P,2);
     272      371902 :   if (l == 3) return Q;
     273      371902 :   if (h == -1)
     274      318735 :     for (i = 3; i < l; i++)
     275             :     {
     276      308186 :       gel(Q,i) = negi(gel(P,i));
     277      308186 :       if (++i == l) break;
     278      304519 :       gel(Q,i) = gel(P,i);
     279             :     }
     280             :   else
     281             :   {
     282             :     GEN hi;
     283      357686 :     gel(Q,3) = mulis(gel(P,3), h);
     284      357686 :     hi = sqrs(h);
     285     1725234 :     for (i = 4; i < l; i++)
     286             :     {
     287     1367548 :       gel(Q,i) = mulii(gel(P,i), hi);
     288     1367548 :       if (i != l-1) hi = mulis(hi,h);
     289             :     }
     290             :   }
     291      371902 :   return Q;
     292             : }
     293             : /* P a ZX, h a t_INT. Return P(h * x), not memory clean; optimize for h = -1 */
     294             : GEN
     295        7378 : ZX_unscale(GEN P, GEN h)
     296             : {
     297             :   long i, l;
     298             :   GEN Q, hi;
     299        7378 :   i = itos_or_0(h); if (i) return ZX_z_unscale(P, i);
     300          14 :   l = lg(P); Q = cgetg(l, t_POL);
     301          14 :   Q[1] = P[1];
     302          14 :   if (l == 2) return Q;
     303          14 :   gel(Q,2) = gel(P,2);
     304          14 :   if (l == 3) return Q;
     305          14 :   hi = h;
     306          14 :   gel(Q,3) = mulii(gel(P,3), hi);
     307          91 :   for (i = 4; i < l; i++)
     308             :   {
     309          77 :     hi = mulii(hi,h);
     310          77 :     gel(Q,i) = mulii(gel(P,i), hi);
     311             :   }
     312          14 :   return Q;
     313             : }
     314             : /* P a ZX. Return P(x << n), not memory clean */
     315             : GEN
     316       62996 : ZX_unscale2n(GEN P, long n)
     317             : {
     318       62996 :   long i, ni = n, l = lg(P);
     319       62996 :   GEN Q = cgetg(l, t_POL);
     320       62996 :   Q[1] = P[1];
     321       62996 :   if (l == 2) return Q;
     322       62996 :   gel(Q,2) = gel(P,2);
     323       62996 :   if (l == 3) return Q;
     324       62996 :   gel(Q,3) = shifti(gel(P,3), ni);
     325      245597 :   for (i=4; i<l; i++)
     326             :   {
     327      182601 :     ni += n;
     328      182601 :     gel(Q,i) = shifti(gel(P,i), ni);
     329             :   }
     330       62996 :   return Q;
     331             : }
     332             : /* P(h*X) / h, assuming h | P(0), i.e. the result is a ZX */
     333             : GEN
     334        1190 : ZX_unscale_div(GEN P, GEN h)
     335             : {
     336        1190 :   long i, l = lg(P);
     337        1190 :   GEN hi, Q = cgetg(l, t_POL);
     338        1190 :   Q[1] = P[1];
     339        1190 :   if (l == 2) return Q;
     340        1190 :   gel(Q,2) = diviiexact(gel(P,2), h);
     341        1190 :   if (l == 3) return Q;
     342        1190 :   gel(Q,3) = gel(P,3);
     343        1190 :   if (l == 4) return Q;
     344        1190 :   hi = h;
     345        1190 :   gel(Q,4) = mulii(gel(P,4), hi);
     346        5257 :   for (i=5; i<l; i++)
     347             :   {
     348        4067 :     hi = mulii(hi,h);
     349        4067 :     gel(Q,i) = mulii(gel(P,i), hi);
     350             :   }
     351        1190 :   return Q;
     352             : }
     353             : 
     354             : GEN
     355         651 : RgXV_unscale(GEN v, GEN h)
     356             : {
     357             :   long i, l;
     358             :   GEN w;
     359         651 :   if (!h || isint1(h)) return v;
     360         350 :   w = cgetg_copy(v, &l);
     361         350 :   for (i=1; i<l; i++) gel(w,i) = RgX_unscale(gel(v,i), h);
     362         350 :   return w;
     363             : }
     364             : 
     365             : /* Return h^degpol(P) P(x / h), not memory clean */
     366             : GEN
     367        2065 : RgX_rescale(GEN P, GEN h)
     368             : {
     369        2065 :   long i, l = lg(P);
     370        2065 :   GEN Q = cgetg(l,t_POL), hi = h;
     371        2065 :   gel(Q,l-1) = gel(P,l-1);
     372        9387 :   for (i=l-2; i>=2; i--)
     373             :   {
     374        9387 :     gel(Q,i) = gmul(gel(P,i), hi);
     375        9387 :     if (i == 2) break;
     376        7322 :     hi = gmul(hi,h);
     377             :   }
     378        2065 :   Q[1] = P[1]; return Q;
     379             : }
     380             : 
     381             : /* A(X^d) --> A(X) */
     382             : GEN
     383      147605 : RgX_deflate(GEN x0, long d)
     384             : {
     385             :   GEN z, y, x;
     386      147605 :   long i,id, dy, dx = degpol(x0);
     387      147605 :   if (d == 1 || dx <= 0) return leafcopy(x0);
     388       78204 :   dy = dx/d;
     389       78204 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     390       78204 :   z = y + 2;
     391       78204 :   x = x0+ 2;
     392       78204 :   for (i=id=0; i<=dy; i++,id+=d) gel(z,i) = gel(x,id);
     393       78204 :   return y;
     394             : }
     395             : 
     396             : /* F a t_RFRAC */
     397             : long
     398          63 : rfrac_deflate_order(GEN F)
     399             : {
     400          63 :   GEN N = gel(F,1), D = gel(F,2);
     401          63 :   long m = (degpol(D) <= 0)? 0: RgX_deflate_order(D);
     402          63 :   if (m == 1) return 1;
     403          28 :   if (typ(N) == t_POL && varn(N) == varn(D))
     404          14 :     m = cgcd(m, RgX_deflate_order(N));
     405          28 :   return m;
     406             : }
     407             : /* F a t_RFRAC */
     408             : GEN
     409          63 : rfrac_deflate_max(GEN F, long *m)
     410             : {
     411          63 :   *m = rfrac_deflate_order(F);
     412          63 :   return rfrac_deflate(F, *m);
     413             : }
     414             : /* F a t_RFRAC */
     415             : GEN
     416          63 : rfrac_deflate(GEN F, long m)
     417             : {
     418          63 :   GEN N = gel(F,1), D = gel(F,2);
     419          63 :   if (m == 1) return F;
     420          28 :   if (typ(N) == t_POL && varn(N) == varn(D)) N = RgX_deflate(N, m);
     421          28 :   D = RgX_deflate(D, m); return mkrfrac(N, D);
     422             : }
     423             : 
     424             : 
     425             : /* return x0(X^d) */
     426             : GEN
     427      593435 : RgX_inflate(GEN x0, long d)
     428             : {
     429      593435 :   long i, id, dy, dx = degpol(x0);
     430      593435 :   GEN x = x0 + 2, z, y;
     431      593435 :   if (dx <= 0) return leafcopy(x0);
     432      586176 :   dy = dx*d;
     433      586176 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     434      586177 :   z = y + 2;
     435      586177 :   for (i=0; i<=dy; i++) gel(z,i) = gen_0;
     436      586177 :   for (i=id=0; i<=dx; i++,id+=d) gel(z,id) = gel(x,i);
     437      586177 :   return y;
     438             : }
     439             : 
     440             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     441             : GEN
     442     1414987 : RgX_translate(GEN P, GEN c)
     443             : {
     444     1414987 :   pari_sp av = avma;
     445             :   GEN Q, *R;
     446             :   long i, k, n;
     447             : 
     448     1414987 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     449     1412534 :   Q = leafcopy(P);
     450     1412534 :   R = (GEN*)(Q+2); n = degpol(P);
     451     1412534 :   if (equali1(c))
     452             :   {
     453        9023 :     for (i=1; i<=n; i++)
     454             :     {
     455        6657 :       for (k=n-i; k<n; k++) R[k] = gadd(R[k], R[k+1]);
     456        6657 :       if (gc_needed(av,2))
     457             :       {
     458           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL(1), i = %ld/%ld", i,n);
     459           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     460             :       }
     461             :     }
     462             :   }
     463     1410168 :   else if (equalim1(c))
     464             :   {
     465      150231 :     for (i=1; i<=n; i++)
     466             :     {
     467      126943 :       for (k=n-i; k<n; k++) R[k] = gsub(R[k], R[k+1]);
     468      126943 :       if (gc_needed(av,2))
     469             :       {
     470           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL(-1), i = %ld/%ld", i,n);
     471           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     472             :       }
     473             :     }
     474             :   }
     475             :   else
     476             :   {
     477     4767176 :     for (i=1; i<=n; i++)
     478             :     {
     479     3380296 :       for (k=n-i; k<n; k++) R[k] = gadd(R[k], gmul(c, R[k+1]));
     480     3380296 :       if (gc_needed(av,2))
     481             :       {
     482           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL, i = %ld/%ld", i,n);
     483           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     484             :       }
     485             :     }
     486             :   }
     487     1412534 :   return gerepilecopy(av, Q);
     488             : }
     489             : 
     490             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     491             : GEN
     492      429425 : ZX_translate(GEN P, GEN c)
     493             : {
     494      429425 :   pari_sp av = avma;
     495             :   GEN Q, *R;
     496             :   long i, k, n;
     497             : 
     498      429425 :   if (!signe(P) || !signe(c)) return ZX_copy(P);
     499      429390 :   Q = leafcopy(P);
     500      429390 :   R = (GEN*)(Q+2); n = degpol(P);
     501      429390 :   if (equali1(c))
     502             :   {
     503     2687982 :     for (i=1; i<=n; i++)
     504             :     {
     505     2363843 :       for (k=n-i; k<n; k++) R[k] = addii(R[k], R[k+1]);
     506     2363843 :       if (gc_needed(av,2))
     507             :       {
     508           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(1), i = %ld/%ld", i,n);
     509           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     510             :       }
     511             :     }
     512             :   }
     513      105251 :   else if (equalim1(c))
     514             :   {
     515          70 :     for (i=1; i<=n; i++)
     516             :     {
     517          49 :       for (k=n-i; k<n; k++) R[k] = subii(R[k], R[k+1]);
     518          49 :       if (gc_needed(av,2))
     519             :       {
     520           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(-1), i = %ld/%ld", i,n);
     521           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     522             :       }
     523             :     }
     524             :   }
     525             :   else
     526             :   {
     527      785792 :     for (i=1; i<=n; i++)
     528             :     {
     529      680562 :       for (k=n-i; k<n; k++) R[k] = addmulii_inplace(R[k], c, R[k+1]);
     530      680562 :       if (gc_needed(av,2))
     531             :       {
     532           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate, i = %ld/%ld", i,n);
     533           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     534             :       }
     535             :     }
     536             :   }
     537      429390 :   return gerepilecopy(av, Q);
     538             : }
     539             : /* return lift( P(X + c) ) using Horner, c in R[y]/(T) */
     540             : GEN
     541       21798 : RgXQX_translate(GEN P, GEN c, GEN T)
     542             : {
     543       21798 :   pari_sp av = avma;
     544             :   GEN Q, *R;
     545             :   long i, k, n;
     546             : 
     547       21798 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     548       21490 :   Q = leafcopy(P);
     549       21490 :   R = (GEN*)(Q+2); n = degpol(P);
     550       74501 :   for (i=1; i<=n; i++)
     551             :   {
     552      222565 :     for (k=n-i; k<n; k++)
     553             :     {
     554      169554 :       pari_sp av2 = avma;
     555      169554 :       R[k] = gerepileupto(av2, RgX_rem(gadd(R[k], gmul(c, R[k+1])), T));
     556             :     }
     557       53011 :     if (gc_needed(av,2))
     558             :     {
     559           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXQX_translate, i = %ld/%ld", i,n);
     560           0 :       Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     561             :     }
     562             :   }
     563       21490 :   return gerepilecopy(av, Q);
     564             : }
     565             : 
     566             : /********************************************************************/
     567             : /**                                                                **/
     568             : /**                          CONVERSIONS                           **/
     569             : /**                       (not memory clean)                       **/
     570             : /**                                                                **/
     571             : /********************************************************************/
     572             : /* to INT / FRAC / (POLMOD mod T), not memory clean because T not copied,
     573             :  * but everything else is */
     574             : static GEN
     575       25106 : QXQ_to_mod_copy(GEN x, GEN T)
     576             : {
     577             :   long d;
     578       25106 :   switch(typ(x))
     579             :   {
     580       10311 :     case t_INT:  return icopy(x);
     581         805 :     case t_FRAC: return gcopy(x);
     582             :     case t_POL:
     583       13990 :       d = degpol(x);
     584       13990 :       if (d < 0) return gen_0;
     585       13430 :       if (d == 0) return gcopy(gel(x,2));
     586       13080 :       return mkpolmod(RgX_copy(x), T);
     587           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     588             :              return NULL;/* LCOV_EXCL_LINE */
     589             :   }
     590             : }
     591             : /* pure shallow version */
     592             : static GEN
     593     1152840 : QXQ_to_mod(GEN x, GEN T)
     594             : {
     595             :   long d;
     596     1152840 :   switch(typ(x))
     597             :   {
     598             :     case t_INT:
     599      438690 :     case t_FRAC: return x;
     600             :     case t_POL:
     601      714150 :       d = degpol(x);
     602      714150 :       if (d < 0) return gen_0;
     603      538996 :       if (d == 0) return gel(x,2);
     604      518018 :       return mkpolmod(x, T);
     605           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     606             :              return NULL;/* LCOV_EXCL_LINE */
     607             :   }
     608             : }
     609             : /* T a ZX, z lifted from (Q[Y]/(T(Y)))[X], apply QXQ_to_mod_copy to all coeffs.
     610             :  * Not memory clean because T not copied, but everything else is */
     611             : static GEN
     612        3416 : QXQX_to_mod(GEN z, GEN T)
     613             : {
     614        3416 :   long i,l = lg(z);
     615        3416 :   GEN x = cgetg(l,t_POL);
     616        3416 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod_copy(gel(z,i), T);
     617        3416 :   x[1] = z[1]; return normalizepol_lg(x,l);
     618             : }
     619             : /* pure shallow version */
     620             : GEN
     621      152711 : QXQX_to_mod_shallow(GEN z, GEN T)
     622             : {
     623      152711 :   long i,l = lg(z);
     624      152711 :   GEN x = cgetg(l,t_POL);
     625      152711 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod(gel(z,i), T);
     626      152711 :   x[1] = z[1]; return normalizepol_lg(x,l);
     627             : }
     628             : /* Apply QXQX_to_mod to all entries. Memory-clean ! */
     629             : GEN
     630        1372 : QXQXV_to_mod(GEN V, GEN T)
     631             : {
     632        1372 :   long i, l = lg(V);
     633        1372 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     634        1372 :   for (i=1;i<l; i++) gel(z,i) = QXQX_to_mod(gel(V,i), T);
     635        1372 :   return z;
     636             : }
     637             : /* Apply QXQ_to_mod_copy to all entries. Memory-clean ! */
     638             : GEN
     639        5181 : QXQV_to_mod(GEN V, GEN T)
     640             : {
     641        5181 :   long i, l = lg(V);
     642        5181 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     643        5181 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod_copy(gel(V,i), T);
     644        5181 :   return z;
     645             : }
     646             : 
     647             : /* Apply QXQ_to_mod_copy to all entries. Memory-clean ! */
     648             : GEN
     649       23926 : QXQC_to_mod_shallow(GEN V, GEN T)
     650             : {
     651       23926 :   long i, l = lg(V);
     652       23926 :   GEN z = cgetg(l, t_COL);
     653       23926 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod(gel(V,i), T);
     654       23926 :   return z;
     655             : }
     656             : 
     657             : GEN
     658        7476 : QXQM_to_mod_shallow(GEN V, GEN T)
     659             : {
     660        7476 :   long i, l = lg(V);
     661        7476 :   GEN z = cgetg(l, t_MAT);
     662        7476 :   for (i=1; i<l; i++) gel(z,i) = QXQC_to_mod_shallow(gel(V,i), T);
     663        7476 :   return z;
     664             : }
     665             : 
     666             : GEN
     667      796615 : RgX_renormalize_lg(GEN x, long lx)
     668             : {
     669             :   long i;
     670     2162197 :   for (i = lx-1; i>1; i--)
     671     2018667 :     if (! gequal0(gel(x,i))) break; /* _not_ isexactzero */
     672      796615 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
     673      796615 :   setlg(x, i+1); setsigne(x, i != 1); return x;
     674             : }
     675             : 
     676             : GEN
     677      521142 : RgV_to_RgX(GEN x, long v)
     678             : {
     679      521142 :   long i, k = lg(x);
     680             :   GEN p;
     681             : 
     682      521142 :   while (--k && gequal0(gel(x,k)));
     683      521142 :   if (!k) return pol_0(v);
     684      515192 :   i = k+2; p = cgetg(i,t_POL);
     685      515192 :   p[1] = evalsigne(1) | evalvarn(v);
     686      515192 :   x--; for (k=2; k<i; k++) gel(p,k) = gel(x,k);
     687      515192 :   return p;
     688             : }
     689             : GEN
     690      158603 : RgV_to_RgX_reverse(GEN x, long v)
     691             : {
     692      158603 :   long j, k, l = lg(x);
     693             :   GEN p;
     694             : 
     695      158603 :   for (k = 1; k < l; k++)
     696      158603 :     if (!gequal0(gel(x,k))) break;
     697      158603 :   if (k == l) return pol_0(v);
     698      158603 :   k -= 1;
     699      158603 :   l -= k;
     700      158603 :   x += k;
     701      158603 :   p = cgetg(l+1,t_POL);
     702      158603 :   p[1] = evalsigne(1) | evalvarn(v);
     703      158603 :   for (j=2, k=l; j<=l; j++) gel(p,j) = gel(x,--k);
     704      158603 :   return p;
     705             : }
     706             : 
     707             : /* return the (N-dimensional) vector of coeffs of p */
     708             : GEN
     709     7107988 : RgX_to_RgC(GEN x, long N)
     710             : {
     711             :   long i, l;
     712             :   GEN z;
     713     7107988 :   l = lg(x)-1; x++;
     714     7107988 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
     715     7107988 :   z = cgetg(N+1,t_COL);
     716     7108061 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
     717     7108061 :   for (   ; i<=N; i++) gel(z,i) = gen_0;
     718     7108061 :   return z;
     719             : }
     720             : GEN
     721      759913 : Rg_to_RgC(GEN x, long N)
     722             : {
     723      759913 :   return (typ(x) == t_POL)? RgX_to_RgC(x,N): scalarcol_shallow(x, N);
     724             : }
     725             : 
     726             : /* vector of polynomials (in v) whose coeffs are given by the columns of x */
     727             : GEN
     728       58112 : RgM_to_RgXV(GEN x, long v)
     729       58112 : { pari_APPLY_type(t_VEC, RgV_to_RgX(gel(x,i), v)) }
     730             : 
     731             : /* matrix whose entries are given by the coeffs of the polynomials in
     732             :  * vector v (considered as degree n-1 polynomials) */
     733             : GEN
     734      135384 : RgV_to_RgM(GEN x, long n)
     735      135384 : { pari_APPLY_type(t_MAT, Rg_to_RgC(gel(x,i), n)) }
     736             : 
     737             : GEN
     738        5454 : RgXV_to_RgM(GEN x, long n)
     739        5454 : { pari_APPLY_type(t_MAT, RgX_to_RgC(gel(x,i), n)) }
     740             : 
     741             : /* polynomial (in v) of polynomials (in w) whose coeffs are given by the columns of x */
     742             : GEN
     743       16845 : RgM_to_RgXX(GEN x, long v,long w)
     744             : {
     745       16845 :   long j, lx = lg(x);
     746       16845 :   GEN y = cgetg(lx+1, t_POL);
     747       16845 :   y[1] = evalsigne(1) | evalvarn(v);
     748       16845 :   y++;
     749       16845 :   for (j=1; j<lx; j++) gel(y,j) = RgV_to_RgX(gel(x,j), w);
     750       16845 :   return normalizepol_lg(--y, lx+1);
     751             : }
     752             : 
     753             : /* matrix whose entries are given by the coeffs of the polynomial v in
     754             :  * two variables (considered as degree n-1 polynomials) */
     755             : GEN
     756         308 : RgXX_to_RgM(GEN v, long n)
     757             : {
     758         308 :   long j, N = lg(v)-1;
     759         308 :   GEN y = cgetg(N, t_MAT);
     760         308 :   for (j=1; j<N; j++) gel(y,j) = Rg_to_RgC(gel(v,j+1), n);
     761         308 :   return y;
     762             : }
     763             : 
     764             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     765             : GEN
     766       19572 : RgXY_swapspec(GEN x, long n, long w, long nx)
     767             : {
     768       19572 :   long j, ly = n+3;
     769       19572 :   GEN y = cgetg(ly, t_POL);
     770       19572 :   y[1] = evalsigne(1);
     771      248426 :   for (j=2; j<ly; j++)
     772             :   {
     773             :     long k;
     774      228854 :     GEN a = cgetg(nx+2,t_POL);
     775      228854 :     a[1] = evalsigne(1) | evalvarn(w);
     776     1169632 :     for (k=0; k<nx; k++)
     777             :     {
     778      940778 :       GEN xk = gel(x,k);
     779      940778 :       if (typ(xk)==t_POL)
     780      865521 :         gel(a,k+2) = j<lg(xk)? gel(xk,j): gen_0;
     781             :       else
     782       75257 :         gel(a,k+2) = j==2 ? xk: gen_0;
     783             :     }
     784      228854 :     gel(y,j) = normalizepol_lg(a, nx+2);
     785             :   }
     786       19572 :   return normalizepol_lg(y,ly);
     787             : }
     788             : 
     789             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     790             : GEN
     791         952 : RgXY_swap(GEN x, long n, long w)
     792             : {
     793         952 :   GEN z = RgXY_swapspec(x+2, n, w, lgpol(x));
     794         952 :   setvarn(z, varn(x)); return z;
     795             : }
     796             : 
     797             : long
     798         176 : RgXY_degreex(GEN b)
     799             : {
     800         176 :   long deg = -1, i;
     801         176 :   if (!signe(b)) return -1;
     802         990 :   for (i = 2; i < lg(b); ++i)
     803             :   {
     804         814 :     GEN bi = gel(b, i);
     805         814 :     if (typ(bi) == t_POL)
     806         785 :       deg = maxss(deg, degpol(bi));
     807             :   }
     808         176 :   return deg;
     809             : }
     810             : 
     811             : /* return (x % X^n). Shallow */
     812             : GEN
     813     2252863 : RgXn_red_shallow(GEN a, long n)
     814             : {
     815     2252863 :   long i, L = n+2, l = lg(a);
     816             :   GEN  b;
     817     2252863 :   if (L >= l) return a; /* deg(x) < n */
     818     1554466 :   b = cgetg(L, t_POL); b[1] = a[1];
     819     1554466 :   for (i=2; i<L; i++) gel(b,i) = gel(a,i);
     820     1554466 :   return normalizepol_lg(b,L);
     821             : }
     822             : 
     823             : GEN
     824         357 : RgXnV_red_shallow(GEN x, long n)
     825         357 : { pari_APPLY_type(t_VEC, RgXn_red_shallow(gel(x,i), n)) }
     826             : 
     827             : /* return (x * X^n). Shallow */
     828             : GEN
     829    64325254 : RgX_shift_shallow(GEN a, long n)
     830             : {
     831    64325254 :   long i, l = lg(a);
     832             :   GEN  b;
     833    64325254 :   if (l == 2 || !n) return a;
     834    48681562 :   l += n;
     835    48681562 :   if (n < 0)
     836             :   {
     837    42091378 :     if (l <= 2) return pol_0(varn(a));
     838    42067746 :     b = cgetg(l, t_POL); b[1] = a[1];
     839    42067748 :     a -= n;
     840    42067748 :     for (i=2; i<l; i++) gel(b,i) = gel(a,i);
     841             :   } else {
     842     6590184 :     b = cgetg(l, t_POL); b[1] = a[1];
     843     6590185 :     a -= n; n += 2;
     844     6590185 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     845     6590185 :     for (   ; i<l; i++) gel(b,i) = gel(a,i);
     846             :   }
     847    48657933 :   return b;
     848             : }
     849             : /* return (x * X^n). */
     850             : GEN
     851     7433141 : RgX_shift(GEN a, long n)
     852             : {
     853     7433141 :   long i, l = lg(a);
     854             :   GEN  b;
     855     7433141 :   if (l == 2 || !n) return RgX_copy(a);
     856     7432875 :   l += n;
     857     7432875 :   if (n < 0)
     858             :   {
     859         840 :     if (l <= 2) return pol_0(varn(a));
     860         798 :     b = cgetg(l, t_POL); b[1] = a[1];
     861         798 :     a -= n;
     862         798 :     for (i=2; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     863             :   } else {
     864     7432035 :     b = cgetg(l, t_POL); b[1] = a[1];
     865     7432035 :     a -= n; n += 2;
     866     7432035 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     867     7432035 :     for (   ; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     868             :   }
     869     7432833 :   return b;
     870             : }
     871             : 
     872             : GEN
     873      316113 : RgX_rotate_shallow(GEN P, long k, long p)
     874             : {
     875      316113 :   long i, l = lgpol(P);
     876             :   GEN r;
     877      316113 :   if (signe(P)==0)
     878         623 :     return pol_0(varn(P));
     879      315490 :   r = cgetg(p+2,t_POL); r[1] = P[1];
     880     2094904 :   for(i=0; i<p; i++)
     881             :   {
     882     1779414 :     long s = 2+(i+k)%p;
     883     1779414 :     gel(r,s) = i<l? gel(P,2+i): gen_0;
     884             :   }
     885      315490 :   return RgX_renormalize(r);
     886             : }
     887             : 
     888             : GEN
     889     2974947 : RgX_mulXn(GEN x, long d)
     890             : {
     891             :   pari_sp av;
     892             :   GEN z;
     893             :   long v;
     894     2974947 :   if (d >= 0) return RgX_shift(x, d);
     895     1463394 :   d = -d;
     896     1463394 :   v = RgX_val(x);
     897     1463394 :   if (v >= d) return RgX_shift(x, -d);
     898     1463380 :   av = avma;
     899     1463380 :   z = gred_rfrac_simple(RgX_shift_shallow(x, -v), pol_xn(d - v, varn(x)));
     900     1463380 :   return gerepileupto(av, z);
     901             : }
     902             : 
     903             : long
     904     2489800 : RgX_val(GEN x)
     905             : {
     906     2489800 :   long i, lx = lg(x);
     907     2489800 :   if (lx == 2) return LONG_MAX;
     908     2522546 :   for (i = 2; i < lx; i++)
     909     2522504 :     if (!isexactzero(gel(x,i))) break;
     910     2489660 :   if (i == lx) return LONG_MAX;/* possible with non-rational zeros */
     911     2489618 :   return i - 2;
     912             : }
     913             : long
     914    38422203 : RgX_valrem(GEN x, GEN *Z)
     915             : {
     916    38422203 :   long v, i, lx = lg(x);
     917    38422203 :   if (lx == 2) { *Z = pol_0(varn(x)); return LONG_MAX; }
     918    72988914 :   for (i = 2; i < lx; i++)
     919    72988914 :     if (!isexactzero(gel(x,i))) break;
     920             :   /* possible with non-rational zeros */
     921    38422203 :   if (i == lx) { *Z = pol_0(varn(x)); return LONG_MAX; }
     922    38422203 :   v = i - 2;
     923    38422203 :   *Z = RgX_shift_shallow(x, -v);
     924    38422203 :   return v;
     925             : }
     926             : long
     927       14956 : RgX_valrem_inexact(GEN x, GEN *Z)
     928             : {
     929             :   long v;
     930       14956 :   if (!signe(x)) { if (Z) *Z = pol_0(varn(x)); return LONG_MAX; }
     931       15432 :   for (v = 0;; v++)
     932       15915 :     if (!gequal0(gel(x,2+v))) break;
     933       14949 :   if (Z) *Z = RgX_shift_shallow(x, -v);
     934       14949 :   return v;
     935             : }
     936             : 
     937             : GEN
     938      109186 : RgXQC_red(GEN x, GEN T)
     939      109186 : { pari_APPLY_type(t_COL, grem(gel(x,i), T)) }
     940             : 
     941             : GEN
     942        1106 : RgXQV_red(GEN x, GEN T)
     943        1106 : { pari_APPLY_type(t_VEC, grem(gel(x,i), T)) }
     944             : 
     945             : GEN
     946       14826 : RgXQM_red(GEN x, GEN T)
     947       14826 : { pari_APPLY_same(RgXQC_red(gel(x,i), T)) }
     948             : 
     949             : GEN
     950           0 : RgXQM_mul(GEN P, GEN Q, GEN T)
     951             : {
     952           0 :   return RgXQM_red(RgM_mul(P, Q), T);
     953             : }
     954             : 
     955             : GEN
     956      422019 : RgXQX_red(GEN P, GEN T)
     957             : {
     958      422019 :   long i, l = lg(P);
     959      422019 :   GEN Q = cgetg(l, t_POL);
     960      422019 :   Q[1] = P[1];
     961      422019 :   for (i=2; i<l; i++) gel(Q,i) = grem(gel(P,i), T);
     962      422019 :   return normalizepol_lg(Q, l);
     963             : }
     964             : 
     965             : GEN
     966      247810 : RgX_deriv(GEN x)
     967             : {
     968      247810 :   long i,lx = lg(x)-1;
     969             :   GEN y;
     970             : 
     971      247810 :   if (lx<3) return pol_0(varn(x));
     972      219712 :   y = cgetg(lx,t_POL); gel(y,2) = gcopy(gel(x,3));
     973      219712 :   for (i=3; i<lx ; i++) gel(y,i) = gmulsg(i-1,gel(x,i+1));
     974      219712 :   y[1] = x[1]; return normalizepol_lg(y,i);
     975             : }
     976             : 
     977             : GEN
     978      530852 : RgX_recipspec_shallow(GEN x, long l, long n)
     979             : {
     980             :   long i;
     981      530852 :   GEN z = cgetg(n+2,t_POL);
     982      530856 :   z[1] = 0; z += 2;
     983    17475554 :   for(i=0; i<l; i++)
     984    16944698 :     gel(z,n-i-1) = gel(x,i);
     985      658205 :   for(   ; i<n; i++)
     986      127349 :     gel(z, n-i-1) = gen_0;
     987      530856 :   return normalizepol_lg(z-2,n+2);
     988             : }
     989             : 
     990             : GEN
     991       45326 : RgXn_recip_shallow(GEN P, long n)
     992             : {
     993       45326 :   GEN Q = RgX_recipspec_shallow(P+2, lgpol(P), n);
     994       45326 :   setvarn(Q, varn(P));
     995       45326 :   return Q;
     996             : }
     997             : 
     998             : /* return coefficients s.t x = x_0 X^n + ... + x_n */
     999             : GEN
    1000        2170 : RgX_recip(GEN x)
    1001             : {
    1002             :   long lx, i, j;
    1003        2170 :   GEN y = cgetg_copy(x, &lx);
    1004        2170 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gcopy(gel(x,j));
    1005        2170 :   return normalizepol_lg(y,lx);
    1006             : }
    1007             : /* shallow version */
    1008             : GEN
    1009      433678 : RgX_recip_shallow(GEN x)
    1010             : {
    1011             :   long lx, i, j;
    1012      433678 :   GEN y = cgetg_copy(x, &lx);
    1013      433727 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
    1014      433727 :   return y;
    1015             : }
    1016             : /*******************************************************************/
    1017             : /*                                                                 */
    1018             : /*                      ADDITION / SUBTRACTION                     */
    1019             : /*                                                                 */
    1020             : /*******************************************************************/
    1021             : /* same variable */
    1022             : GEN
    1023    20463687 : RgX_add(GEN x, GEN y)
    1024             : {
    1025    20463687 :   long i, lx = lg(x), ly = lg(y);
    1026             :   GEN z;
    1027    20463687 :   if (ly <= lx) {
    1028    18119270 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1029    18119270 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1030    18119270 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1031    18119270 :     z = normalizepol_lg(z, lx);
    1032             :   } else {
    1033     2344417 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1034     2344417 :     for (i=2; i < lx; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1035     2344417 :     for (   ; i < ly; i++) gel(z,i) = gcopy(gel(y,i));
    1036     2344417 :     z = normalizepol_lg(z, ly);
    1037             :   }
    1038    20463687 :   return z;
    1039             : }
    1040             : GEN
    1041    22287461 : RgX_sub(GEN x, GEN y)
    1042             : {
    1043    22287461 :   long i, lx = lg(x), ly = lg(y);
    1044             :   GEN z;
    1045    22287461 :   if (ly <= lx) {
    1046    18852093 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1047    18852093 :     for (i=2; i < ly; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1048    18852093 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1049    18852093 :     z = normalizepol_lg(z, lx);
    1050             :   } else {
    1051     3435368 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1052     3435368 :     for (i=2; i < lx; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1053     3435368 :     for (   ; i < ly; i++) gel(z,i) = gneg(gel(y,i));
    1054     3435368 :     z = normalizepol_lg(z, ly);
    1055             :   }
    1056    22287461 :   return z;
    1057             : }
    1058             : GEN
    1059     1240502 : RgX_neg(GEN x)
    1060             : {
    1061     1240502 :   long i, lx = lg(x);
    1062     1240502 :   GEN y = cgetg(lx, t_POL); y[1] = x[1];
    1063     1240502 :   for (i=2; i<lx; i++) gel(y,i) = gneg(gel(x,i));
    1064     1240502 :   return y;
    1065             : }
    1066             : 
    1067             : GEN
    1068    12592260 : RgX_Rg_add(GEN y, GEN x)
    1069             : {
    1070             :   GEN z;
    1071    12592260 :   long lz = lg(y), i;
    1072    12592260 :   if (lz == 2) return scalarpol(x,varn(y));
    1073     9610164 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1074     9610166 :   gel(z,2) = gadd(gel(y,2),x);
    1075     9610095 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1076             :   /* probably useless unless lz = 3, but cannot be skipped if y is
    1077             :    * an inexact 0 */
    1078     9610174 :   return normalizepol_lg(z,lz);
    1079             : }
    1080             : GEN
    1081        2830 : RgX_Rg_add_shallow(GEN y, GEN x)
    1082             : {
    1083             :   GEN z;
    1084        2830 :   long lz = lg(y), i;
    1085        2830 :   if (lz == 2) return scalarpol(x,varn(y));
    1086        2830 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1087        2830 :   gel(z,2) = gadd(gel(y,2),x);
    1088        2830 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
    1089        2830 :   return z = normalizepol_lg(z,lz);
    1090             : }
    1091             : GEN
    1092       32699 : RgX_Rg_sub(GEN y, GEN x)
    1093             : {
    1094             :   GEN z;
    1095       32699 :   long lz = lg(y), i;
    1096       32699 :   if (lz == 2)
    1097             :   { /* scalarpol(gneg(x),varn(y)) optimized */
    1098        2576 :     long v = varn(y);
    1099        2576 :     if (isrationalzero(x)) return pol_0(v);
    1100          14 :     z = cgetg(3,t_POL);
    1101          28 :     z[1] = gequal0(x)? evalvarn(v)
    1102          14 :                    : evalvarn(v) | evalsigne(1);
    1103          14 :     gel(z,2) = gneg(x); return z;
    1104             :   }
    1105       30123 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1106       30123 :   gel(z,2) = gsub(gel(y,2),x);
    1107       30123 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1108       30123 :   return z = normalizepol_lg(z,lz);
    1109             : }
    1110             : GEN
    1111      528961 : Rg_RgX_sub(GEN x, GEN y)
    1112             : {
    1113             :   GEN z;
    1114      528961 :   long lz = lg(y), i;
    1115      528961 :   if (lz == 2) return scalarpol(x,varn(y));
    1116      527855 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1117      527855 :   gel(z,2) = gsub(x, gel(y,2));
    1118      527855 :   for(i=3; i<lz; i++) gel(z,i) = gneg(gel(y,i));
    1119      527855 :   return z = normalizepol_lg(z,lz);
    1120             : }
    1121             : /*******************************************************************/
    1122             : /*                                                                 */
    1123             : /*                  KARATSUBA MULTIPLICATION                       */
    1124             : /*                                                                 */
    1125             : /*******************************************************************/
    1126             : #if 0
    1127             : /* to debug Karatsuba-like routines */
    1128             : GEN
    1129             : zx_debug_spec(GEN x, long nx)
    1130             : {
    1131             :   GEN z = cgetg(nx+2,t_POL);
    1132             :   long i;
    1133             :   for (i=0; i<nx; i++) gel(z,i+2) = stoi(x[i]);
    1134             :   z[1] = evalsigne(1); return z;
    1135             : }
    1136             : 
    1137             : GEN
    1138             : RgX_debug_spec(GEN x, long nx)
    1139             : {
    1140             :   GEN z = cgetg(nx+2,t_POL);
    1141             :   long i;
    1142             :   for (i=0; i<nx; i++) z[i+2] = x[i];
    1143             :   z[1] = evalsigne(1); return z;
    1144             : }
    1145             : #endif
    1146             : 
    1147             : /* generic multiplication */
    1148             : GEN
    1149     3387352 : RgX_addspec_shallow(GEN x, GEN y, long nx, long ny)
    1150             : {
    1151             :   GEN z, t;
    1152             :   long i;
    1153     3387352 :   if (nx == ny) {
    1154      772008 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1155      772008 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1156      772008 :     return normalizepol_lg(z, nx+2);
    1157             :   }
    1158     2615344 :   if (ny < nx) {
    1159     2460412 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1160     2460412 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1161     2460412 :     for (   ; i < nx; i++) gel(t,i) = gel(x,i);
    1162     2460412 :     return normalizepol_lg(z, nx+2);
    1163             :   } else {
    1164      154932 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1165      154932 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1166      154932 :     for (   ; i < ny; i++) gel(t,i) = gel(y,i);
    1167      154932 :     return normalizepol_lg(z, ny+2);
    1168             :   }
    1169             : }
    1170             : GEN
    1171      209438 : RgX_addspec(GEN x, GEN y, long nx, long ny)
    1172             : {
    1173             :   GEN z, t;
    1174             :   long i;
    1175      209438 :   if (nx == ny) {
    1176       12425 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1177       12425 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1178       12425 :     return normalizepol_lg(z, nx+2);
    1179             :   }
    1180      197013 :   if (ny < nx) {
    1181      195312 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1182      195312 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1183      195312 :     for (   ; i < nx; i++) gel(t,i) = gcopy(gel(x,i));
    1184      195312 :     return normalizepol_lg(z, nx+2);
    1185             :   } else {
    1186        1701 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1187        1701 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1188        1701 :     for (   ; i < ny; i++) gel(t,i) = gcopy(gel(y,i));
    1189        1701 :     return normalizepol_lg(z, ny+2);
    1190             :   }
    1191             : }
    1192             : 
    1193             : /* Return the vector of coefficients of x, where we replace rational 0s by NULL
    1194             :  * [ to speed up basic operation s += x[i]*y[j] ]. We create a proper
    1195             :  * t_VECSMALL, to hold this, which can be left on stack: gerepile
    1196             :  * will not crash on it. The returned vector itself is not a proper GEN,
    1197             :  * we access the coefficients as x[i], i = 0..deg(x) */
    1198             : static GEN
    1199    21251004 : RgXspec_kill0(GEN x, long lx)
    1200             : {
    1201    21251004 :   GEN z = cgetg(lx+1, t_VECSMALL) + 1; /* inhibit gerepile-wise */
    1202             :   long i;
    1203    86271638 :   for (i=0; i <lx; i++)
    1204             :   {
    1205    65020634 :     GEN c = gel(x,i);
    1206    65020634 :     z[i] = (long)(isrationalzero(c)? NULL: c);
    1207             :   }
    1208    21251004 :   return z;
    1209             : }
    1210             : 
    1211             : INLINE GEN
    1212    47700942 : RgX_mulspec_basecase_limb(GEN x, GEN y, long a, long b)
    1213             : {
    1214    47700942 :   pari_sp av = avma;
    1215    47700942 :   GEN s = NULL;
    1216             :   long i;
    1217             : 
    1218   226990675 :   for (i=a; i<b; i++)
    1219   179289733 :     if (gel(y,i) && gel(x,-i))
    1220             :     {
    1221   101278840 :       GEN t = gmul(gel(y,i), gel(x,-i));
    1222   101278840 :       s = s? gadd(s, t): t;
    1223             :     }
    1224    47700942 :   return s? gerepileupto(av, s): gen_0;
    1225             : }
    1226             : 
    1227             : /* assume nx >= ny > 0, return x * y * t^v */
    1228             : static GEN
    1229     8458671 : RgX_mulspec_basecase(GEN x, GEN y, long nx, long ny, long v)
    1230             : {
    1231             :   long i, lz, nz;
    1232             :   GEN z;
    1233             : 
    1234     8458671 :   x = RgXspec_kill0(x,nx);
    1235     8458671 :   y = RgXspec_kill0(y,ny);
    1236     8458671 :   lz = nx + ny + 1; nz = lz-2;
    1237     8458671 :   lz += v;
    1238     8458671 :   z = cgetg(lz, t_POL) + 2; /* x:y:z [i] = term of degree i */
    1239     8458671 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1240     8458671 :   for (i=0; i<ny; i++)gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0, i+1);
    1241     8458671 :   for (  ; i<nx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ny);
    1242     8458671 :   for (  ; i<nz; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-nx+1,ny);
    1243     8458671 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1244             : }
    1245             : 
    1246             : /* return (x * X^d) + y. Assume d > 0 */
    1247             : GEN
    1248     2521290 : RgX_addmulXn_shallow(GEN x0, GEN y0, long d)
    1249             : {
    1250             :   GEN x, y, xd, yd, zd;
    1251             :   long a, lz, nx, ny;
    1252             : 
    1253     2521290 :   if (!signe(x0)) return y0;
    1254     2495241 :   ny = lgpol(y0);
    1255     2495241 :   nx = lgpol(x0);
    1256     2495241 :   zd = (GEN)avma;
    1257     2495241 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1258     2495241 :   if (a <= 0)
    1259             :   {
    1260      202480 :     lz = nx+d+2;
    1261      202480 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1262      202480 :     while (xd > x) gel(--zd,0) = gel(--xd,0);
    1263      202480 :     x = zd + a;
    1264      202480 :     while (zd > x) gel(--zd,0) = gen_0;
    1265             :   }
    1266             :   else
    1267             :   {
    1268     2292761 :     xd = new_chunk(d); yd = y+d;
    1269     2292761 :     x = RgX_addspec_shallow(x,yd, nx,a);
    1270     2292761 :     lz = (a>nx)? ny+2: lg(x)+d;
    1271     2292761 :     x += 2; while (xd > x) *--zd = *--xd;
    1272             :   }
    1273     2495241 :   while (yd > y) *--zd = *--yd;
    1274     2495241 :   *--zd = x0[1];
    1275     2495241 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1276             : }
    1277             : GEN
    1278      501398 : RgX_addmulXn(GEN x0, GEN y0, long d)
    1279             : {
    1280             :   GEN x, y, xd, yd, zd;
    1281             :   long a, lz, nx, ny;
    1282             : 
    1283      501398 :   if (!signe(x0)) return RgX_copy(y0);
    1284      500821 :   nx = lgpol(x0);
    1285      500821 :   ny = lgpol(y0);
    1286      500821 :   zd = (GEN)avma;
    1287      500821 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1288      500821 :   if (a <= 0)
    1289             :   {
    1290      291383 :     lz = nx+d+2;
    1291      291383 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1292      291383 :     while (xd > x) gel(--zd,0) = gcopy(gel(--xd,0));
    1293      291383 :     x = zd + a;
    1294      291383 :     while (zd > x) gel(--zd,0) = gen_0;
    1295             :   }
    1296             :   else
    1297             :   {
    1298      209438 :     xd = new_chunk(d); yd = y+d;
    1299      209438 :     x = RgX_addspec(x,yd, nx,a);
    1300      209438 :     lz = (a>nx)? ny+2: lg(x)+d;
    1301      209438 :     x += 2; while (xd > x) *--zd = *--xd;
    1302             :   }
    1303      500821 :   while (yd > y) gel(--zd,0) = gcopy(gel(--yd,0));
    1304      500821 :   *--zd = x0[1];
    1305      500821 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1306             : }
    1307             : 
    1308             : /* return x * y mod t^n */
    1309             : static GEN
    1310     2163539 : RgXn_mul_basecase(GEN x, GEN y, long n)
    1311             : {
    1312     2163539 :   long i, lz = n+2, lx = lgpol(x), ly = lgpol(y);
    1313             :   GEN z;
    1314     2163539 :   if (lx < 0) return pol_0(varn(x));
    1315     2163539 :   if (ly < 0) return pol_0(varn(x));
    1316     2163539 :   z = cgetg(lz, t_POL) + 2;
    1317     2163539 :   x+=2; if (lx > n) lx = n;
    1318     2163539 :   y+=2; if (ly > n) ly = n;
    1319     2163539 :   z[-1] = x[-1];
    1320     2163539 :   if (ly > lx) { swap(x,y); lswap(lx,ly); }
    1321     2163539 :   x = RgXspec_kill0(x, lx);
    1322     2163539 :   y = RgXspec_kill0(y, ly);
    1323             :   /* x:y:z [i] = term of degree i */
    1324     2163539 :   for (i=0;i<ly; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,i+1);
    1325     2163539 :   for (  ; i<lx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ly);
    1326     2163539 :   for (  ; i<n; i++)  gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-lx+1,ly);
    1327     2163539 :   return normalizepol_lg(z - 2, lz);
    1328             : }
    1329             : /* Mulders / Karatsuba product f*g mod t^n (Hanrot-Zimmermann variant) */
    1330             : static GEN
    1331     2254098 : RgXn_mul2(GEN f, GEN g, long n)
    1332             : {
    1333     2254098 :   pari_sp av = avma;
    1334             :   GEN fe,fo, ge,go, l,h,m;
    1335             :   long n0, n1;
    1336     2254098 :   if (degpol(f) + degpol(g) < n) return RgX_mul(f,g);
    1337     2193485 :   if (n < 80) return RgXn_mul_basecase(f,g,n);
    1338       29946 :   n0 = n>>1; n1 = n-n0;
    1339       29946 :   RgX_even_odd(f, &fe, &fo);
    1340       29946 :   RgX_even_odd(g, &ge, &go);
    1341       29946 :   l = RgXn_mul(fe,ge,n1);
    1342       29946 :   h = RgXn_mul(fo,go,n0);
    1343       29946 :   m = RgX_sub(RgXn_mul(RgX_add(fe,fo),RgX_add(ge,go),n0), RgX_add(l,h));
    1344             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1345             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1346       29946 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1347             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1348       29946 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1349       29946 :   m = RgX_inflate(m,2);
    1350             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1351       29946 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1352       29946 :   h = RgX_inflate(h,2);
    1353       29946 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1354       29946 :   return gerepileupto(av, h);
    1355             : }
    1356             : /* (f*g) \/ x^n */
    1357             : static GEN
    1358        2974 : RgX_mulhigh_i2(GEN f, GEN g, long n)
    1359             : {
    1360        2974 :   long d = degpol(f)+degpol(g) + 1 - n;
    1361             :   GEN h;
    1362        2974 :   if (d <= 2) return RgX_shift_shallow(RgX_mul(f,g), -n);
    1363        1886 :   h = RgX_recip_shallow(RgXn_mul(RgX_recip_shallow(f),
    1364             :                                  RgX_recip_shallow(g), d));
    1365        1886 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1366             : }
    1367             : 
    1368             : /* (f*g) \/ x^n */
    1369             : static GEN
    1370           0 : RgX_sqrhigh_i2(GEN f, long n)
    1371             : {
    1372           0 :   long d = 2*degpol(f)+ 1 - n;
    1373             :   GEN h;
    1374           0 :   if (d <= 2) return RgX_shift_shallow(RgX_sqr(f), -n);
    1375           0 :   h = RgX_recip_shallow(RgXn_sqr(RgX_recip_shallow(f), d));
    1376           0 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1377             : }
    1378             : 
    1379             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
    1380             :  * b+2 were sent instead. na, nb = number of terms of a, b.
    1381             :  * Only c, c0, c1, c2 are genuine GEN.
    1382             :  */
    1383             : GEN
    1384     9031063 : RgX_mulspec(GEN a, GEN b, long na, long nb)
    1385             : {
    1386             :   GEN a0, c, c0;
    1387     9031063 :   long n0, n0a, i, v = 0;
    1388             :   pari_sp av;
    1389             : 
    1390     9031063 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v++; }
    1391     9031063 :   while (nb && isrationalzero(gel(b,0))) { b++; nb--; v++; }
    1392     9031063 :   if (na < nb) swapspec(a,b, na,nb);
    1393     9031063 :   if (!nb) return pol_0(0);
    1394             : 
    1395     8930084 :   if (nb < RgX_MUL_LIMIT) return RgX_mulspec_basecase(a,b,na,nb, v);
    1396      471413 :   RgX_shift_inplace_init(v);
    1397      471413 :   i = (na>>1); n0 = na-i; na = i;
    1398      471413 :   av = avma; a0 = a+n0; n0a = n0;
    1399      471413 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1400             : 
    1401      471413 :   if (nb > n0)
    1402             :   {
    1403             :     GEN b0,c1,c2;
    1404             :     long n0b;
    1405             : 
    1406      470040 :     nb -= n0; b0 = b+n0; n0b = n0;
    1407      470040 :     while (n0b && isrationalzero(gel(b,n0b-1))) n0b--;
    1408      470040 :     c = RgX_mulspec(a,b,n0a,n0b);
    1409      470040 :     c0 = RgX_mulspec(a0,b0, na,nb);
    1410             : 
    1411      470040 :     c2 = RgX_addspec_shallow(a0,a, na,n0a);
    1412      470040 :     c1 = RgX_addspec_shallow(b0,b, nb,n0b);
    1413             : 
    1414      470040 :     c1 = RgX_mulspec(c1+2,c2+2, lgpol(c1),lgpol(c2));
    1415      470040 :     c2 = RgX_sub(c1, RgX_add(c0,c));
    1416      470040 :     c0 = RgX_addmulXn_shallow(c0, c2, n0);
    1417             :   }
    1418             :   else
    1419             :   {
    1420        1373 :     c = RgX_mulspec(a,b,n0a,nb);
    1421        1373 :     c0 = RgX_mulspec(a0,b,na,nb);
    1422             :   }
    1423      471413 :   c0 = RgX_addmulXn(c0,c,n0);
    1424      471413 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1425             : }
    1426             : 
    1427             : INLINE GEN
    1428       34814 : RgX_sqrspec_basecase_limb(GEN x, long a, long i)
    1429             : {
    1430       34814 :   pari_sp av = avma;
    1431       34814 :   GEN s = NULL;
    1432       34814 :   long j, l = (i+1)>>1;
    1433      101003 :   for (j=a; j<l; j++)
    1434             :   {
    1435       66189 :     GEN xj = gel(x,j), xx = gel(x,i-j);
    1436       66189 :     if (xj && xx)
    1437             :     {
    1438       58867 :       GEN t = gmul(xj, xx);
    1439       58867 :       s = s? gadd(s, t): t;
    1440             :     }
    1441             :   }
    1442       34814 :   if (s) s = gshift(s,1);
    1443       34814 :   if ((i&1) == 0)
    1444             :   {
    1445       20699 :     GEN t = gel(x, i>>1);
    1446       20699 :     if (t) {
    1447       18382 :       t = gsqr(t);
    1448       18382 :       s = s? gadd(s, t): t;
    1449             :     }
    1450             :   }
    1451       34814 :   return s? gerepileupto(av,s): gen_0;
    1452             : }
    1453             : static GEN
    1454        6584 : RgX_sqrspec_basecase(GEN x, long nx, long v)
    1455             : {
    1456             :   long i, lz, nz;
    1457             :   GEN z;
    1458             : 
    1459        6584 :   if (!nx) return pol_0(0);
    1460        6584 :   x = RgXspec_kill0(x,nx);
    1461        6584 :   lz = (nx << 1) + 1, nz = lz-2;
    1462        6584 :   lz += v;
    1463        6584 :   z = cgetg(lz,t_POL) + 2;
    1464        6584 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1465        6584 :   for (i=0; i<nx; i++)gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1466        6584 :   for (  ; i<nz; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, i-nx+1, i);
    1467        6584 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1468             : }
    1469             : /* return x^2 mod t^n */
    1470             : static GEN
    1471           0 : RgXn_sqr_basecase(GEN x, long n)
    1472             : {
    1473           0 :   long i, lz = n+2, lx = lgpol(x);
    1474             :   GEN z;
    1475           0 :   if (lx < 0) return pol_0(varn(x));
    1476           0 :   z = cgetg(lz, t_POL);
    1477           0 :   z[1] = x[1];
    1478           0 :   x+=2; if (lx > n) lx = n;
    1479           0 :   x = RgXspec_kill0(x,lx);
    1480           0 :   z+=2;/* x:z [i] = term of degree i */
    1481           0 :   for (i=0;i<lx; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1482           0 :   for (  ; i<n; i++)  gel(z,i) = RgX_sqrspec_basecase_limb(x, i-lx+1, i);
    1483           0 :   z -= 2; return normalizepol_lg(z, lz);
    1484             : }
    1485             : /* Mulders / Karatsuba product f^2 mod t^n (Hanrot-Zimmermann variant) */
    1486             : static GEN
    1487         238 : RgXn_sqr2(GEN f, long n)
    1488             : {
    1489         238 :   pari_sp av = avma;
    1490             :   GEN fe,fo, l,h,m;
    1491             :   long n0, n1;
    1492         238 :   if (2*degpol(f) < n) return RgX_sqr_i(f);
    1493           0 :   if (n < 80) return RgXn_sqr_basecase(f,n);
    1494           0 :   n0 = n>>1; n1 = n-n0;
    1495           0 :   RgX_even_odd(f, &fe, &fo);
    1496           0 :   l = RgXn_sqr(fe,n1);
    1497           0 :   h = RgXn_sqr(fo,n0);
    1498           0 :   m = RgX_sub(RgXn_sqr(RgX_add(fe,fo),n0), RgX_add(l,h));
    1499             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1500             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1501           0 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1502             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1503           0 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1504           0 :   m = RgX_inflate(m,2);
    1505             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1506           0 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1507           0 :   h = RgX_inflate(h,2);
    1508           0 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1509           0 :   return gerepileupto(av, h);
    1510             : }
    1511             : GEN
    1512        6623 : RgX_sqrspec(GEN a, long na)
    1513             : {
    1514             :   GEN a0, c, c0, c1;
    1515        6623 :   long n0, n0a, i, v = 0;
    1516             :   pari_sp av;
    1517             : 
    1518        6623 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v += 2; }
    1519        6623 :   if (na<RgX_SQR_LIMIT) return RgX_sqrspec_basecase(a, na, v);
    1520          39 :   RgX_shift_inplace_init(v);
    1521          39 :   i = (na>>1); n0 = na-i; na = i;
    1522          39 :   av = avma; a0 = a+n0; n0a = n0;
    1523          39 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1524             : 
    1525          39 :   c = RgX_sqrspec(a,n0a);
    1526          39 :   c0 = RgX_sqrspec(a0,na);
    1527          39 :   c1 = gmul2n(RgX_mulspec(a0,a, na,n0a), 1);
    1528          39 :   c0 = RgX_addmulXn_shallow(c0,c1, n0);
    1529          39 :   c0 = RgX_addmulXn(c0,c,n0);
    1530          39 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1531             : }
    1532             : 
    1533             : /* (X^a + A)(X^b + B) - X^(a+b), where deg A < a, deg B < b */
    1534             : GEN
    1535      748919 : RgX_mul_normalized(GEN A, long a, GEN B, long b)
    1536             : {
    1537      748919 :   GEN z = RgX_mul(A, B);
    1538      748919 :   if (a < b)
    1539        5446 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(A, B, b-a), z, a);
    1540      743473 :   else if (a > b)
    1541      449852 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(B, A, a-b), z, b);
    1542             :   else
    1543      293621 :     z = RgX_addmulXn_shallow(RgX_add(A, B), z, a);
    1544      748919 :   return z;
    1545             : }
    1546             : 
    1547             : GEN
    1548     7618158 : RgX_mul_i(GEN x, GEN y)
    1549             : {
    1550     7618158 :   GEN z = RgX_mulspec(x+2, y+2, lgpol(x), lgpol(y));
    1551     7618158 :   setvarn(z, varn(x)); return z;
    1552             : }
    1553             : 
    1554             : GEN
    1555        6545 : RgX_sqr_i(GEN x)
    1556             : {
    1557        6545 :   GEN z = RgX_sqrspec(x+2, lgpol(x));
    1558        6545 :   setvarn(z,varn(x)); return z;
    1559             : }
    1560             : 
    1561             : /*******************************************************************/
    1562             : /*                                                                 */
    1563             : /*                               DIVISION                          */
    1564             : /*                                                                 */
    1565             : /*******************************************************************/
    1566             : GEN
    1567      456833 : RgX_Rg_divexact(GEN x, GEN y) {
    1568             :   long i, lx;
    1569             :   GEN z;
    1570      456833 :   if (typ(y) == t_INT && is_pm1(y))
    1571       31857 :     return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1572      424976 :   z = cgetg_copy(x, &lx); z[1] = x[1];
    1573      424976 :   for (i=2; i<lx; i++) gel(z,i) = gdivexact(gel(x,i),y);
    1574      424976 :   return z;
    1575             : }
    1576             : GEN
    1577    23949179 : RgX_Rg_div(GEN x, GEN y) {
    1578             :   long i, lx;
    1579    23949179 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1580    23949179 :   for (i=2; i<lx; i++) gel(z,i) = gdiv(gel(x,i),y);
    1581    23949179 :   return normalizepol_lg(z, lx);
    1582             : }
    1583             : GEN
    1584       13223 : RgX_normalize(GEN x)
    1585             : {
    1586       13223 :   GEN z, d = NULL;
    1587       13223 :   long i, n = lg(x)-1;
    1588       13223 :   for (i = n; i > 1; i--) { d = gel(x,i); if (!gequal0(d)) break; }
    1589       13223 :   if (i == 1) return pol_0(varn(x));
    1590       13223 :   if (i == n && isint1(d)) return x;
    1591        6020 :   n = i; z = cgetg(n+1, t_POL); z[1] = x[1];
    1592        6020 :   for (i=2; i<n; i++) gel(z,i) = gdiv(gel(x,i),d);
    1593        6020 :   gel(z,n) = Rg_get_1(d); return z;
    1594             : }
    1595             : GEN
    1596        2625 : RgX_divs(GEN x, long y) {
    1597             :   long i, lx;
    1598        2625 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1599        2625 :   for (i=2; i<lx; i++) gel(z,i) = gdivgs(gel(x,i),y);
    1600        2625 :   return normalizepol_lg(z, lx);
    1601             : }
    1602             : GEN
    1603       50634 : RgX_div_by_X_x(GEN a, GEN x, GEN *r)
    1604             : {
    1605       50634 :   long l = lg(a), i;
    1606       50634 :   GEN a0, z0, z = cgetg(l-1, t_POL);
    1607       50634 :   z[1] = a[1];
    1608       50634 :   a0 = a + l-1;
    1609       50634 :   z0 = z + l-2; *z0 = *a0--;
    1610      871445 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1611             :   {
    1612      820811 :     GEN t = gadd(gel(a0--,0), gmul(x, gel(z0--,0)));
    1613      820811 :     gel(z0,0) = t;
    1614             :   }
    1615       50634 :   if (r) *r = gadd(gel(a0,0), gmul(x, gel(z0,0)));
    1616       50634 :   return z;
    1617             : }
    1618             : /* Polynomial division x / y:
    1619             :  *   if pr = ONLY_REM return remainder, otherwise return quotient
    1620             :  *   if pr = ONLY_DIVIDES return quotient if division is exact, else NULL
    1621             :  *   if pr != NULL set *pr to remainder, as the last object on stack */
    1622             : /* assume, typ(x) = typ(y) = t_POL, same variable */
    1623             : static GEN
    1624    12445190 : RgX_divrem_i(GEN x, GEN y, GEN *pr)
    1625             : {
    1626             :   pari_sp avy, av, av1;
    1627             :   long dx,dy,dz,i,j,sx,lr;
    1628             :   GEN z,p1,p2,rem,y_lead,mod,p;
    1629             :   GEN (*f)(GEN,GEN);
    1630             : 
    1631    12445190 :   if (!signe(y)) pari_err_INV("RgX_divrem",y);
    1632             : 
    1633    12445190 :   dy = degpol(y);
    1634    12445190 :   y_lead = gel(y,dy+2);
    1635    12445190 :   if (gequal0(y_lead)) /* normalize denominator if leading term is 0 */
    1636             :   {
    1637           0 :     pari_warn(warner,"normalizing a polynomial with 0 leading term");
    1638           0 :     for (dy--; dy>=0; dy--)
    1639             :     {
    1640           0 :       y_lead = gel(y,dy+2);
    1641           0 :       if (!gequal0(y_lead)) break;
    1642             :     }
    1643             :   }
    1644    12445190 :   if (!dy) /* y is constant */
    1645             :   {
    1646        3016 :     if (pr == ONLY_REM) return pol_0(varn(x));
    1647        3016 :     z = RgX_Rg_div(x, y_lead);
    1648        3016 :     if (pr == ONLY_DIVIDES) return z;
    1649        1196 :     if (pr) *pr = pol_0(varn(x));
    1650        1196 :     return z;
    1651             :   }
    1652    12442174 :   dx = degpol(x);
    1653    12442174 :   if (dx < dy)
    1654             :   {
    1655      662965 :     if (pr == ONLY_REM) return RgX_copy(x);
    1656      371733 :     if (pr == ONLY_DIVIDES) return signe(x)? NULL: pol_0(varn(x));
    1657      371712 :     z = pol_0(varn(x));
    1658      371712 :     if (pr) *pr = RgX_copy(x);
    1659      371712 :     return z;
    1660             :   }
    1661             : 
    1662             :   /* x,y in R[X], y non constant */
    1663    11779209 :   av = avma;
    1664    11779209 :   p = NULL;
    1665    11779209 :   if (RgX_is_FpX(x, &p) && RgX_is_FpX(y, &p) && p)
    1666             :   {
    1667      225841 :     z = FpX_divrem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p, pr);
    1668      225841 :     if (!z) return gc_NULL(av);
    1669      225841 :     z = FpX_to_mod(z, p);
    1670      225841 :     if (!pr || pr == ONLY_REM || pr == ONLY_DIVIDES)
    1671      121681 :       return gerepileupto(av, z);
    1672      104160 :     *pr = FpX_to_mod(*pr, p);
    1673      104160 :     gerepileall(av, 2, pr, &z);
    1674      104160 :     return z;
    1675             :   }
    1676    11553368 :   switch(typ(y_lead))
    1677             :   {
    1678             :     case t_REAL:
    1679           0 :       y_lead = ginv(y_lead);
    1680           0 :       f = gmul; mod = NULL;
    1681           0 :       break;
    1682             :     case t_INTMOD:
    1683        2246 :     case t_POLMOD: y_lead = ginv(y_lead);
    1684        2246 :       f = gmul; mod = gmodulo(gen_1, gel(y_lead,1));
    1685        2246 :       break;
    1686    11551122 :     default: if (gequal1(y_lead)) y_lead = NULL;
    1687    11551122 :       f = gdiv; mod = NULL;
    1688             :   }
    1689             : 
    1690    11553368 :   if (y_lead == NULL)
    1691     9768759 :     p2 = gel(x,dx+2);
    1692             :   else {
    1693             :     for(;;) {
    1694     1784609 :       p2 = f(gel(x,dx+2),y_lead);
    1695     1784609 :       p2 = simplify_shallow(p2);
    1696     1784609 :       if (!isexactzero(p2) || (--dx < 0)) break;
    1697             :     }
    1698     1784609 :     if (dx < dy) /* leading coeff of x was in fact zero */
    1699             :     {
    1700           0 :       if (pr == ONLY_DIVIDES) {
    1701           0 :         set_avma(av);
    1702           0 :         return (dx < 0)? pol_0(varn(x)) : NULL;
    1703             :       }
    1704           0 :       if (pr == ONLY_REM)
    1705             :       {
    1706           0 :         if (dx < 0)
    1707           0 :           return gerepilecopy(av, scalarpol(p2, varn(x)));
    1708             :         else
    1709             :         {
    1710             :           GEN t;
    1711           0 :           set_avma(av);
    1712           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1713           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1714           0 :           return t;
    1715             :         }
    1716             :       }
    1717           0 :       if (pr) /* cf ONLY_REM above */
    1718             :       {
    1719           0 :         if (dx < 0)
    1720             :         {
    1721           0 :           p2 = gclone(p2);
    1722           0 :           set_avma(av);
    1723           0 :           z = pol_0(varn(x));
    1724           0 :           x = scalarpol(p2, varn(x));
    1725           0 :           gunclone(p2);
    1726             :         }
    1727             :         else
    1728             :         {
    1729             :           GEN t;
    1730           0 :           set_avma(av);
    1731           0 :           z = pol_0(varn(x));
    1732           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1733           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1734           0 :           x = t;
    1735             :         }
    1736           0 :         *pr = x;
    1737             :       }
    1738             :       else
    1739             :       {
    1740           0 :         set_avma(av);
    1741           0 :         z = pol_0(varn(x));
    1742             :       }
    1743           0 :       return z;
    1744             :     }
    1745             :   }
    1746             :   /* dx >= dy */
    1747    11553368 :   avy = avma;
    1748    11553368 :   dz = dx-dy;
    1749    11553368 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1750    11553368 :   x += 2;
    1751    11553368 :   z += 2;
    1752    11553368 :   y += 2;
    1753    11553368 :   gel(z,dz) = gcopy(p2);
    1754             : 
    1755    28822242 :   for (i=dx-1; i>=dy; i--)
    1756             :   {
    1757    17268874 :     av1=avma; p1=gel(x,i);
    1758    17268874 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1759    17268874 :     if (y_lead) p1 = simplify(f(p1,y_lead));
    1760             : 
    1761    17268874 :     if (isrationalzero(p1)) { set_avma(av1); p1 = gen_0; }
    1762             :     else
    1763     7717237 :       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
    1764    17268874 :     gel(z,i-dy) = p1;
    1765             :   }
    1766    11553368 :   if (!pr) return gerepileupto(av,z-2);
    1767             : 
    1768     4868966 :   rem = (GEN)avma; av1 = (pari_sp)new_chunk(dx+3);
    1769     5228640 :   for (sx=0; ; i--)
    1770             :   {
    1771     5588314 :     p1 = gel(x,i);
    1772             :     /* we always enter this loop at least once */
    1773     5228640 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1774     5228640 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1775     5228640 :     if (!gequal0(p1)) { sx = 1; break; } /* remainder is non-zero */
    1776     3430244 :     if (!isexactzero(p1)) break;
    1777     3416158 :     if (!i) break;
    1778      359674 :     set_avma(av1);
    1779             :   }
    1780     4868966 :   if (pr == ONLY_DIVIDES)
    1781             :   {
    1782        1659 :     if (sx) return gc_NULL(av);
    1783        1652 :     set_avma((pari_sp)rem); return gerepileupto(av,z-2);
    1784             :   }
    1785     4867307 :   lr=i+3; rem -= lr;
    1786     4867307 :   if (avma==av1) { set_avma((pari_sp)rem); p1 = gcopy(p1); }
    1787     4823695 :   else p1 = gerepileupto((pari_sp)rem,p1);
    1788     4867307 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1789     4867307 :   rem[1] = z[-1];
    1790     4867307 :   rem += 2;
    1791     4867307 :   gel(rem,i) = p1;
    1792     6175184 :   for (i--; i>=0; i--)
    1793             :   {
    1794     1307877 :     av1=avma; p1 = gel(x,i);
    1795     1307877 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1796     1307877 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1797     1307877 :     gel(rem,i) = avma==av1? gcopy(p1):gerepileupto(av1,p1);
    1798             :   }
    1799     4867307 :   rem -= 2;
    1800     4867307 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1801     4867307 :   if (pr == ONLY_REM) return gerepileupto(av,rem);
    1802     4111060 :   z -= 2;
    1803             :   {
    1804     4111060 :     GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
    1805     4111060 :     gerepilemanysp(av,avy,gptr,2); *pr = rem; return z;
    1806             :   }
    1807             : }
    1808             : 
    1809             : GEN
    1810    11397711 : RgX_divrem(GEN x, GEN y, GEN *pr)
    1811             : {
    1812    11397711 :   if (pr == ONLY_REM) return RgX_rem(x, y);
    1813    11397711 :   return RgX_divrem_i(x, y, pr);
    1814             : }
    1815             : 
    1816             : /* x and y in (R[Y]/T)[X]  (lifted), T in R[Y]. y preferably monic */
    1817             : GEN
    1818       74214 : RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)
    1819             : {
    1820             :   long vx, dx, dy, dz, i, j, sx, lr;
    1821             :   pari_sp av0, av, tetpil;
    1822             :   GEN z,p1,rem,lead;
    1823             : 
    1824       74214 :   if (!signe(y)) pari_err_INV("RgXQX_divrem",y);
    1825       74214 :   vx = varn(x);
    1826       74214 :   dx = degpol(x);
    1827       74214 :   dy = degpol(y);
    1828       74214 :   if (dx < dy)
    1829             :   {
    1830       24213 :     if (pr)
    1831             :     {
    1832       24213 :       av0 = avma; x = RgXQX_red(x, T);
    1833       24213 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: gen_0; }
    1834       24213 :       if (pr == ONLY_REM) return x;
    1835           0 :       *pr = x;
    1836             :     }
    1837           0 :     return pol_0(vx);
    1838             :   }
    1839       50001 :   lead = leading_coeff(y);
    1840       50001 :   if (!dy) /* y is constant */
    1841             :   {
    1842         546 :     if (pr && pr != ONLY_DIVIDES)
    1843             :     {
    1844           0 :       if (pr == ONLY_REM) return pol_0(vx);
    1845           0 :       *pr = pol_0(vx);
    1846             :     }
    1847         546 :     if (gequal1(lead)) return RgX_copy(x);
    1848           0 :     av0 = avma; x = gmul(x, ginvmod(lead,T)); tetpil = avma;
    1849           0 :     return gerepile(av0,tetpil,RgXQX_red(x,T));
    1850             :   }
    1851       49455 :   av0 = avma; dz = dx-dy;
    1852       49455 :   lead = gequal1(lead)? NULL: gclone(ginvmod(lead,T));
    1853       49455 :   set_avma(av0);
    1854       49455 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1855       49455 :   x += 2; y += 2; z += 2;
    1856             : 
    1857       49455 :   p1 = gel(x,dx); av = avma;
    1858       49455 :   gel(z,dz) = lead? gerepileupto(av, grem(gmul(p1,lead), T)): gcopy(p1);
    1859      293205 :   for (i=dx-1; i>=dy; i--)
    1860             :   {
    1861      243750 :     av=avma; p1=gel(x,i);
    1862      243750 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1863      243750 :     if (lead) p1 = gmul(grem(p1, T), lead);
    1864      243750 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil, grem(p1, T));
    1865             :   }
    1866       49455 :   if (!pr) { guncloneNULL(lead); return z-2; }
    1867             : 
    1868       48209 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    1869       78231 :   for (sx=0; ; i--)
    1870             :   {
    1871      108253 :     p1 = gel(x,i);
    1872       78231 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1873       78231 :     tetpil=avma; p1 = grem(p1, T); if (!gequal0(p1)) { sx = 1; break; }
    1874       35630 :     if (!i) break;
    1875       30022 :     set_avma(av);
    1876             :   }
    1877       48209 :   if (pr == ONLY_DIVIDES)
    1878             :   {
    1879        2605 :     guncloneNULL(lead);
    1880        2605 :     if (sx) return gc_NULL(av0);
    1881        2500 :     set_avma((pari_sp)rem); return z-2;
    1882             :   }
    1883       45604 :   lr=i+3; rem -= lr;
    1884       45604 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1885       45604 :   rem[1] = z[-1];
    1886       45604 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    1887       45604 :   rem += 2; gel(rem,i) = p1;
    1888      120230 :   for (i--; i>=0; i--)
    1889             :   {
    1890       74626 :     av=avma; p1 = gel(x,i);
    1891      277026 :     for (j=0; j<=i && j<=dz; j++)
    1892      202400 :       p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1893       74626 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, grem(p1, T));
    1894             :   }
    1895       45604 :   rem -= 2;
    1896       45604 :   guncloneNULL(lead);
    1897       45604 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1898       45604 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    1899         119 :   *pr = rem; return z-2;
    1900             : }
    1901             : 
    1902             : /*******************************************************************/
    1903             : /*                                                                 */
    1904             : /*                        PSEUDO-DIVISION                          */
    1905             : /*                                                                 */
    1906             : /*******************************************************************/
    1907             : INLINE GEN
    1908      523875 : rem(GEN c, GEN T)
    1909             : {
    1910      523875 :   if (T && typ(c) == t_POL && varn(c) == varn(T)) c = RgX_rem(c, T);
    1911      523875 :   return c;
    1912             : }
    1913             : 
    1914             : /* x, y, are ZYX, lc(y) is an integer, T is a ZY */
    1915             : int
    1916        7210 : ZXQX_dvd(GEN x, GEN y, GEN T)
    1917             : {
    1918             :   long dx, dy, dz, i, p, T_ismonic;
    1919        7210 :   pari_sp av = avma, av2;
    1920             :   GEN y_lead;
    1921             : 
    1922        7210 :   if (!signe(y)) pari_err_INV("ZXQX_dvd",y);
    1923        7210 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1924        7210 :   if (typ(y_lead) == t_POL) y_lead = gel(y_lead, 2); /* t_INT */
    1925             :   /* if monic, no point in using pseudo-division */
    1926        7210 :   if (gequal1(y_lead)) return signe(RgXQX_rem(x, y, T)) == 0;
    1927        4207 :   T_ismonic = gequal1(leading_coeff(T));
    1928        4207 :   dx = degpol(x);
    1929        4207 :   if (dx < dy) return !signe(x);
    1930        4207 :   (void)new_chunk(2);
    1931        4207 :   x = RgX_recip_shallow(x)+2;
    1932        4207 :   y = RgX_recip_shallow(y)+2;
    1933             :   /* pay attention to sparse divisors */
    1934        8892 :   for (i = 1; i <= dy; i++)
    1935        4685 :     if (!signe(gel(y,i))) gel(y,i) = NULL;
    1936        4207 :   dz = dx-dy; p = dz+1;
    1937        4207 :   av2 = avma;
    1938             :   for (;;)
    1939       32233 :   {
    1940       36440 :     GEN m, x0 = gel(x,0), y0 = y_lead, cx = content(x0);
    1941       36440 :     x0 = gneg(x0); p--;
    1942       36440 :     m = gcdii(cx, y0);
    1943       36440 :     if (!equali1(m))
    1944             :     {
    1945       35320 :       x0 = gdiv(x0, m);
    1946       35320 :       y0 = diviiexact(y0, m);
    1947       35320 :       if (equali1(y0)) y0 = NULL;
    1948             :     }
    1949       76587 :     for (i=1; i<=dy; i++)
    1950             :     {
    1951       40147 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1952       40147 :       if (gel(y,i)) c = gadd(c, gmul(x0,gel(y,i)));
    1953       40147 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1954       40147 :       gel(x,i) = c;
    1955             :     }
    1956      410083 :     for (   ; i<=dx; i++)
    1957             :     {
    1958      373643 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1959      373643 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1960      373643 :       gel(x,i) = c;
    1961             :     }
    1962       40987 :     do { x++; dx--; } while (dx >= 0 && !signe(gel(x,0)));
    1963       36440 :     if (dx < dy) break;
    1964       32233 :     if (gc_needed(av2,1))
    1965             :     {
    1966           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZXQX_dvd dx = %ld >= %ld",dx,dy);
    1967           0 :       gerepilecoeffs(av2,x,dx+1);
    1968             :     }
    1969             :   }
    1970        4207 :   return gc_bool(av, dx < 0);
    1971             : }
    1972             : 
    1973             : /* T either NULL or a t_POL. */
    1974             : GEN
    1975       29546 : RgXQX_pseudorem(GEN x, GEN y, GEN T)
    1976             : {
    1977       29546 :   long vx = varn(x), dx, dy, dz, i, lx, p;
    1978       29546 :   pari_sp av = avma, av2;
    1979             :   GEN y_lead;
    1980             : 
    1981       29546 :   if (!signe(y)) pari_err_INV("RgXQX_pseudorem",y);
    1982       29546 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1983             :   /* if monic, no point in using pseudo-division */
    1984       29546 :   if (gequal1(y_lead)) return T? RgXQX_rem(x, y, T): RgX_rem(x, y);
    1985       26074 :   dx = degpol(x);
    1986       26074 :   if (dx < dy) return RgX_copy(x);
    1987       26074 :   (void)new_chunk(2);
    1988       26074 :   x = RgX_recip_shallow(x)+2;
    1989       26074 :   y = RgX_recip_shallow(y)+2;
    1990             :   /* pay attention to sparse divisors */
    1991       81195 :   for (i = 1; i <= dy; i++)
    1992       55121 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    1993       26074 :   dz = dx-dy; p = dz+1;
    1994       26074 :   av2 = avma;
    1995             :   for (;;)
    1996             :   {
    1997      158622 :     gel(x,0) = gneg(gel(x,0)); p--;
    1998      295371 :     for (i=1; i<=dy; i++)
    1999             :     {
    2000      203023 :       GEN c = gmul(y_lead, gel(x,i));
    2001      203023 :       if (gel(y,i)) c = gadd(c, gmul(gel(x,0),gel(y,i)));
    2002      203023 :       gel(x,i) = rem(c, T);
    2003             :     }
    2004      313527 :     for (   ; i<=dx; i++)
    2005             :     {
    2006      221179 :       GEN c = gmul(y_lead, gel(x,i));
    2007      221179 :       gel(x,i) = rem(c, T);
    2008             :     }
    2009      101932 :     do { x++; dx--; } while (dx >= 0 && gequal0(gel(x,0)));
    2010       92348 :     if (dx < dy) break;
    2011       66274 :     if (gc_needed(av2,1))
    2012             :     {
    2013           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudorem dx = %ld >= %ld",dx,dy);
    2014           0 :       gerepilecoeffs(av2,x,dx+1);
    2015             :     }
    2016             :   }
    2017       26074 :   if (dx < 0) return pol_0(vx);
    2018       24135 :   lx = dx+3; x -= 2;
    2019       24135 :   x[0] = evaltyp(t_POL) | evallg(lx);
    2020       24135 :   x[1] = evalsigne(1) | evalvarn(vx);
    2021       24135 :   x = RgX_recip_shallow(x);
    2022       24135 :   if (p)
    2023             :   { /* multiply by y[0]^p   [beware dummy vars from FpX_FpXY_resultant] */
    2024        2668 :     GEN t = y_lead;
    2025        2668 :     if (T && typ(t) == t_POL && varn(t) == varn(T))
    2026           0 :       t = RgXQ_powu(t, p, T);
    2027             :     else
    2028        2668 :       t = gpowgs(t, p);
    2029        7447 :     for (i=2; i<lx; i++)
    2030             :     {
    2031        4779 :       GEN c = gmul(gel(x,i), t);
    2032        4779 :       gel(x,i) = rem(c,T);
    2033             :     }
    2034        2668 :     if (!T) return gerepileupto(av, x);
    2035             :   }
    2036       21467 :   return gerepilecopy(av, x);
    2037             : }
    2038             : 
    2039             : GEN
    2040       29546 : RgX_pseudorem(GEN x, GEN y) { return RgXQX_pseudorem(x,y, NULL); }
    2041             : 
    2042             : /* Compute z,r s.t lc(y)^(dx-dy+1) x = z y + r */
    2043             : GEN
    2044        7769 : RgXQX_pseudodivrem(GEN x, GEN y, GEN T, GEN *ptr)
    2045             : {
    2046        7769 :   long vx = varn(x), dx, dy, dz, i, iz, lx, lz, p;
    2047        7769 :   pari_sp av = avma, av2;
    2048             :   GEN z, r, ypow, y_lead;
    2049             : 
    2050        7769 :   if (!signe(y)) pari_err_INV("RgXQX_pseudodivrem",y);
    2051        7769 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2052        7769 :   if (gequal1(y_lead)) return T? RgXQX_divrem(x,y, T, ptr): RgX_divrem(x,y, ptr);
    2053        7223 :   dx = degpol(x);
    2054        7223 :   if (dx < dy) { *ptr = RgX_copy(x); return pol_0(vx); }
    2055        7223 :   if (dx == dy)
    2056             :   {
    2057          70 :     GEN x_lead = gel(x,lg(x)-1);
    2058          70 :     x = RgX_renormalize_lg(leafcopy(x), lg(x)-1);
    2059          70 :     y = RgX_renormalize_lg(leafcopy(y), lg(y)-1);
    2060          70 :     r = RgX_sub(RgX_Rg_mul(x, y_lead), RgX_Rg_mul(y, x_lead));
    2061          70 :     *ptr = gerepileupto(av, r); return scalarpol(x_lead, vx);
    2062             :   }
    2063        7153 :   (void)new_chunk(2);
    2064        7153 :   x = RgX_recip_shallow(x)+2;
    2065        7153 :   y = RgX_recip_shallow(y)+2;
    2066             :   /* pay attention to sparse divisors */
    2067       36628 :   for (i = 1; i <= dy; i++)
    2068       29475 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2069        7153 :   dz = dx-dy; p = dz+1;
    2070        7153 :   lz = dz+3;
    2071        7153 :   z = cgetg(lz, t_POL);
    2072        7153 :   z[1] = evalsigne(1) | evalvarn(vx);
    2073        7153 :   for (i = 2; i < lz; i++) gel(z,i) = gen_0;
    2074        7153 :   ypow = new_chunk(dz+1);
    2075        7153 :   gel(ypow,0) = gen_1;
    2076        7153 :   gel(ypow,1) = y_lead;
    2077       13131 :   for (i=2; i<=dz; i++)
    2078             :   {
    2079        5978 :     GEN c = gmul(gel(ypow,i-1), y_lead);
    2080        5978 :     gel(ypow,i) = rem(c,T);
    2081             :   }
    2082        7153 :   av2 = avma;
    2083        7153 :   for (iz=2;;)
    2084             :   {
    2085       22743 :     p--;
    2086       14948 :     gel(z,iz++) = rem(gmul(gel(x,0), gel(ypow,p)), T);
    2087       64831 :     for (i=1; i<=dy; i++)
    2088             :     {
    2089       49883 :       GEN c = gmul(y_lead, gel(x,i));
    2090       49883 :       if (gel(y,i)) c = gsub(c, gmul(gel(x,0),gel(y,i)));
    2091       49883 :       gel(x,i) = rem(c, T);
    2092             :     }
    2093       39033 :     for (   ; i<=dx; i++)
    2094             :     {
    2095       24085 :       GEN c = gmul(y_lead, gel(x,i));
    2096       24085 :       gel(x,i) = rem(c,T);
    2097             :     }
    2098       14948 :     x++; dx--;
    2099       14948 :     while (dx >= dy && gequal0(gel(x,0))) { x++; dx--; iz++; }
    2100       14948 :     if (dx < dy) break;
    2101        7795 :     if (gc_needed(av2,1))
    2102             :     {
    2103           0 :       GEN X = x-2;
    2104           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudodivrem dx=%ld >= %ld",dx,dy);
    2105           0 :       X[0] = evaltyp(t_POL)|evallg(dx+3); X[1] = z[1]; /* hack */
    2106           0 :       gerepileall(av2,2, &X, &z); x = X+2;
    2107             :     }
    2108             :   }
    2109        7153 :   while (dx >= 0 && gequal0(gel(x,0))) { x++; dx--; }
    2110        7153 :   if (dx < 0)
    2111          98 :     x = pol_0(vx);
    2112             :   else
    2113             :   {
    2114        7055 :     lx = dx+3; x -= 2;
    2115        7055 :     x[0] = evaltyp(t_POL) | evallg(lx);
    2116        7055 :     x[1] = evalsigne(1) | evalvarn(vx);
    2117        7055 :     x = RgX_recip_shallow(x);
    2118             :   }
    2119        7153 :   z = RgX_recip_shallow(z);
    2120        7153 :   r = x;
    2121        7153 :   if (p)
    2122             :   {
    2123        3207 :     GEN c = gel(ypow,p); r = RgX_Rg_mul(r, c);
    2124        3207 :     if (T && typ(c) == t_POL && varn(c) == varn(T)) r = RgXQX_red(r, T);
    2125             :   }
    2126        7153 :   gerepileall(av, 2, &z, &r);
    2127        7153 :   *ptr = r; return z;
    2128             : }
    2129             : GEN
    2130        7573 : RgX_pseudodivrem(GEN x, GEN y, GEN *ptr)
    2131        7573 : { return RgXQX_pseudodivrem(x,y,NULL,ptr); }
    2132             : 
    2133             : GEN
    2134       12789 : RgXQX_mul(GEN x, GEN y, GEN T)
    2135             : {
    2136       12789 :   return RgXQX_red(RgX_mul(x,y), T);
    2137             : }
    2138             : GEN
    2139    70325211 : RgX_Rg_mul(GEN y, GEN x) {
    2140             :   long i, ly;
    2141    70325211 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2142    70325211 :   if (ly == 2) return z;
    2143    70239762 :   for (i = 2; i < ly; i++) gel(z,i) = gmul(x,gel(y,i));
    2144    70239755 :   return normalizepol_lg(z,ly);
    2145             : }
    2146             : GEN
    2147       14700 : RgX_muls(GEN y, long x) {
    2148             :   long i, ly;
    2149       14700 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2150       14700 :   if (ly == 2) return z;
    2151       14602 :   for (i = 2; i < ly; i++) gel(z,i) = gmulsg(x,gel(y,i));
    2152       14602 :   return normalizepol_lg(z,ly);
    2153             : }
    2154             : GEN
    2155        4753 : RgXQX_RgXQ_mul(GEN x, GEN y, GEN T)
    2156             : {
    2157        4753 :   return RgXQX_red(RgX_Rg_mul(x,y), T);
    2158             : }
    2159             : GEN
    2160          77 : RgXQV_RgXQ_mul(GEN v, GEN x, GEN T)
    2161             : {
    2162          77 :   return RgXQV_red(RgV_Rg_mul(v,x), T);
    2163             : }
    2164             : 
    2165             : GEN
    2166        1792 : RgXQX_sqr(GEN x, GEN T)
    2167             : {
    2168        1792 :   return RgXQX_red(RgX_sqr(x), T);
    2169             : }
    2170             : 
    2171             : GEN
    2172         448 : RgXQX_powers(GEN P, long n, GEN T)
    2173             : {
    2174         448 :   GEN v = cgetg(n+2, t_VEC);
    2175             :   long i;
    2176         448 :   gel(v, 1) = pol_1(varn(T));
    2177         448 :   if (n==0) return v;
    2178         448 :   gel(v, 2) = gcopy(P);
    2179         448 :   for (i = 2; i <= n; i++) gel(v,i+1) = RgXQX_mul(P, gel(v,i), T);
    2180         448 :   return v;
    2181             : }
    2182             : 
    2183             : static GEN
    2184      144171 : _add(void *data, GEN x, GEN y) { (void)data; return RgX_add(x, y); }
    2185             : static GEN
    2186           0 : _sub(void *data, GEN x, GEN y) { (void)data; return RgX_sub(x, y); }
    2187             : static GEN
    2188       56300 : _sqr(void *data, GEN x) { return RgXQ_sqr(x, (GEN)data); }
    2189             : static GEN
    2190       46474 : _mul(void *data, GEN x, GEN y) { return RgXQ_mul(x,y, (GEN)data); }
    2191             : static GEN
    2192      261777 : _cmul(void *data, GEN P, long a, GEN x) { (void)data; return RgX_Rg_mul(x,gel(P,a+2)); }
    2193             : static GEN
    2194      259705 : _one(void *data) { return pol_1(varn((GEN)data)); }
    2195             : static GEN
    2196         553 : _zero(void *data) { return pol_0(varn((GEN)data)); }
    2197             : static GEN
    2198      156253 : _red(void *data, GEN x) { (void)data; return gcopy(x); }
    2199             : 
    2200             : static struct bb_algebra RgXQ_algebra = { _red, _add, _sub,
    2201             :               _mul, _sqr, _one, _zero };
    2202             : 
    2203             : GEN
    2204           0 : RgX_RgXQV_eval(GEN Q, GEN x, GEN T)
    2205             : {
    2206           0 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)T,&RgXQ_algebra,_cmul);
    2207             : }
    2208             : 
    2209             : GEN
    2210      117438 : RgX_RgXQ_eval(GEN Q, GEN x, GEN T)
    2211             : {
    2212      117438 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2213      117438 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)T,&RgXQ_algebra,_cmul);
    2214             : }
    2215             : 
    2216             : /* mod X^n */
    2217             : struct modXn {
    2218             :   long v; /* varn(X) */
    2219             :   long n;
    2220             : } ;
    2221             : static GEN
    2222        2268 : _sqrXn(void *data, GEN x) {
    2223        2268 :   struct modXn *S = (struct modXn*)data;
    2224        2268 :   return RgXn_sqr(x, S->n);
    2225             : }
    2226             : static GEN
    2227        1309 : _mulXn(void *data, GEN x, GEN y) {
    2228        1309 :   struct modXn *S = (struct modXn*)data;
    2229        1309 :   return RgXn_mul(x,y, S->n);
    2230             : }
    2231             : static GEN
    2232        1477 : _oneXn(void *data) {
    2233        1477 :   struct modXn *S = (struct modXn*)data;
    2234        1477 :   return pol_1(S->v);
    2235             : }
    2236             : static GEN
    2237           0 : _zeroXn(void *data) {
    2238           0 :   struct modXn *S = (struct modXn*)data;
    2239           0 :   return pol_0(S->v);
    2240             : }
    2241             : static struct bb_algebra RgXn_algebra = { _red, _add, _sub, _mulXn, _sqrXn,
    2242             :                                           _oneXn, _zeroXn };
    2243             : 
    2244             : GEN
    2245         357 : RgXn_powers(GEN x, long m, long n)
    2246             : {
    2247         357 :   long d = degpol(x);
    2248         357 :   int use_sqr = (d<<1) >= n;
    2249             :   struct modXn S;
    2250         357 :   S.v = varn(x); S.n = n;
    2251         357 :   return gen_powers(x,m,use_sqr,(void*)&S,_sqrXn,_mulXn,_oneXn);
    2252             : }
    2253             : 
    2254             : GEN
    2255        1848 : RgXn_powu_i(GEN x, ulong m, long n)
    2256             : {
    2257             :   struct modXn S;
    2258        1848 :   S.v = varn(x); S.n = n;
    2259        1848 :   return gen_powu_i(x, m, (void*)&S,_sqrXn,_mulXn);
    2260             : }
    2261             : GEN
    2262           0 : RgXn_powu(GEN x, ulong m, long n)
    2263             : {
    2264             :   struct modXn S;
    2265           0 :   S.v = varn(x); S.n = n;
    2266           0 :   return gen_powu(x, m, (void*)&S,_sqrXn,_mulXn);
    2267             : }
    2268             : 
    2269             : GEN
    2270         714 : RgX_RgXnV_eval(GEN Q, GEN x, long n)
    2271             : {
    2272             :   struct modXn S;
    2273         714 :   S.v = varn(gel(x,2)); S.n = n;
    2274         714 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)&S,&RgXn_algebra,_cmul);
    2275             : }
    2276             : 
    2277             : GEN
    2278           0 : RgX_RgXn_eval(GEN Q, GEN x, long n)
    2279             : {
    2280           0 :   int use_sqr = 2*degpol(x) >= n;
    2281             :   struct modXn S;
    2282           0 :   S.v = varn(x); S.n = n;
    2283           0 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2284             : }
    2285             : 
    2286             : /* Q(x) mod t^n, x in R[t], n >= 1 */
    2287             : GEN
    2288        6629 : RgXn_eval(GEN Q, GEN x, long n)
    2289             : {
    2290        6629 :   long d = degpol(x);
    2291             :   int use_sqr;
    2292             :   struct modXn S;
    2293        6629 :   if (d == 1 && isrationalzero(gel(x,2)))
    2294             :   {
    2295        6622 :     GEN y = RgX_unscale(Q, gel(x,3));
    2296        6622 :     setvarn(y, varn(x)); return y;
    2297             :   }
    2298           7 :   S.v = varn(x);
    2299           7 :   S.n = n;
    2300           7 :   use_sqr = (d<<1) >= n;
    2301           7 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2302             : }
    2303             : 
    2304             : /* (f*g mod t^n) \ t^n2, assuming 2*n2 >= n */
    2305             : static GEN
    2306       44813 : RgXn_mulhigh(GEN f, GEN g, long n2, long n)
    2307             : {
    2308       44813 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2309       44813 :   return RgX_add(RgX_mulhigh_i(fl, g, n2), RgXn_mul(fh, g, n - n2));
    2310             : }
    2311             : 
    2312             : /* (f^2 mod t^n) \ t^n2, assuming 2*n2 >= n */
    2313             : static GEN
    2314           0 : RgXn_sqrhigh(GEN f, long n2, long n)
    2315             : {
    2316           0 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2317           0 :   return RgX_add(RgX_mulhigh_i(fl, f, n2), RgXn_mul(fh, f, n - n2));
    2318             : }
    2319             : 
    2320             : GEN
    2321      126516 : RgXn_inv_i(GEN f, long e)
    2322             : {
    2323             :   pari_sp av;
    2324             :   ulong mask;
    2325             :   GEN W, a;
    2326      126516 :   long v = varn(f), n = 1;
    2327             : 
    2328      126516 :   if (!signe(f)) pari_err_INV("RgXn_inv",f);
    2329      126516 :   a = ginv(gel(f,2));
    2330      126508 :   if (e == 1) return scalarpol(a, v);
    2331      126480 :   else if (e == 2)
    2332             :   {
    2333             :     GEN b;
    2334      108812 :     if (degpol(f) <= 0 || gequal0(b = gel(f,3))) return scalarpol(a, v);
    2335       88767 :     av = avma; b = gneg(b);
    2336       88767 :     if (!gequal1(a)) b = gmul(b, gsqr(a));
    2337       88767 :     W = deg1pol_shallow(b, a, v);
    2338       88766 :     return gcopy(W);
    2339             :   }
    2340       17668 :   W = scalarpol_shallow(a,v);
    2341       17668 :   mask = quadratic_prec_mask(e);
    2342       17668 :   av = avma;
    2343       80149 :   for (;mask>1;)
    2344             :   {
    2345             :     GEN u, fr;
    2346       44813 :     long n2 = n;
    2347       44813 :     n<<=1; if (mask & 1) n--;
    2348       44813 :     mask >>= 1;
    2349       44813 :     fr = RgXn_red_shallow(f, n);
    2350       44813 :     u = RgXn_mul(W, RgXn_mulhigh(fr, W, n2, n), n-n2);
    2351       44813 :     W = RgX_sub(W, RgX_shift_shallow(u, n2));
    2352       44813 :     if (gc_needed(av,2))
    2353             :     {
    2354           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_inv, e = %ld", n);
    2355           0 :       W = gerepileupto(av, W);
    2356             :     }
    2357             :   }
    2358       17668 :   return W;
    2359             : }
    2360             : 
    2361             : static GEN
    2362           0 : RgXn_inv_FpX(GEN x, long e, GEN p)
    2363             : {
    2364           0 :   pari_sp av = avma;
    2365             :   GEN r;
    2366           0 :   if (lgefint(p) == 3)
    2367             :   {
    2368           0 :     ulong pp = uel(p, 2);
    2369           0 :     r = Flx_to_ZX_inplace(Flxn_inv(RgX_to_Flx(x, pp), e, pp));
    2370             :   }
    2371             :   else
    2372           0 :     r = FpXn_inv(RgX_to_FpX(x, p), e, p);
    2373           0 :   return gerepileupto(av, FpX_to_mod(r, p));
    2374             : }
    2375             : 
    2376             : static GEN
    2377           0 : RgXn_inv_FpXQX(GEN x, long n, GEN pol, GEN p)
    2378             : {
    2379           0 :   pari_sp av = avma;
    2380           0 :   GEN r, T = RgX_to_FpX(pol, p);
    2381           0 :   if (signe(T) == 0) pari_err_OP("/", gen_1, x);
    2382           0 :   r = FpXQXn_inv(RgX_to_FpXQX(x, T, p), n, T, p);
    2383           0 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2384             : }
    2385             : 
    2386             : #define code(t1,t2) ((t1 << 6) | t2)
    2387             : 
    2388             : static GEN
    2389       83269 : RgXn_inv_fast(GEN x, long e)
    2390             : {
    2391             :   GEN p, pol;
    2392             :   long pa;
    2393       83269 :   long t = RgX_type(x,&p,&pol,&pa);
    2394       83269 :   switch(t)
    2395             :   {
    2396           0 :     case t_INTMOD: return RgXn_inv_FpX(x, e, p);
    2397             :     case code(t_POLMOD, t_INTMOD):
    2398           0 :                    return RgXn_inv_FpXQX(x, e, pol, p);
    2399       83269 :     default:       return NULL;
    2400             :   }
    2401             : }
    2402             : #undef code
    2403             : 
    2404             : GEN
    2405       83270 : RgXn_inv(GEN f, long e)
    2406             : {
    2407       83270 :   GEN h = RgXn_inv_fast(f, e);
    2408       83269 :   if (h) return h;
    2409       83269 :   return RgXn_inv_i(f, e);
    2410             : }
    2411             : 
    2412             : GEN
    2413       12971 : RgXn_exp(GEN h, long e)
    2414             : {
    2415       12971 :   pari_sp av = avma, av2;
    2416       12971 :   long v = varn(h), n=1;
    2417       12971 :   GEN f = pol_1(v), g = pol_1(v);
    2418       12971 :   ulong mask = quadratic_prec_mask(e);
    2419       12971 :   av2 = avma;
    2420       12971 :   if (signe(h)==0 || degpol(h)<1 || !gequal0(gel(h,2)))
    2421           0 :     pari_err_DOMAIN("RgXn_exp","valuation", "<", gen_1, h);
    2422       40152 :   for (;mask>1;)
    2423             :   {
    2424             :     GEN q, w;
    2425       14210 :     long n2 = n;
    2426       14210 :     n<<=1; if (mask & 1) n--;
    2427       14210 :     mask >>= 1;
    2428       14210 :     g = RgX_sub(RgX_muls(g,2),RgXn_mul(f,RgXn_sqr(g,n2),n2));
    2429       14210 :     q = RgX_deriv(RgXn_red_shallow(h,n2));
    2430       14210 :     w = RgX_add(q, RgXn_mul(g, RgX_sub(RgX_deriv(f), RgXn_mul(f,q,n-1)),n-1));
    2431       14210 :     f = RgX_add(f, RgXn_mul(f, RgX_sub(RgXn_red_shallow(h, n), RgX_integ(w)), n));
    2432       14210 :     if (gc_needed(av2,2))
    2433             :     {
    2434           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_exp, e = %ld", n);
    2435           0 :       gerepileall(av2, 2, &f, &g);
    2436             :     }
    2437             :   }
    2438       12971 :   return gerepileupto(av, f);
    2439             : }
    2440             : 
    2441             : GEN
    2442          91 : RgXn_reverse(GEN f, long e)
    2443             : {
    2444          91 :   pari_sp av = avma, av2;
    2445             :   ulong mask;
    2446             :   GEN fi, a, df, W, an;
    2447          91 :   long v = varn(f), n=1;
    2448          91 :   if (degpol(f)<1 || !gequal0(gel(f,2)))
    2449           0 :     pari_err_INV("serreverse",f);
    2450          91 :   fi = ginv(gel(f,3));
    2451          91 :   a = deg1pol_shallow(fi,gen_0,v);
    2452          91 :   if (e <= 2) return gerepilecopy(av, a);
    2453          91 :   W = scalarpol(fi,v);
    2454          91 :   df = RgX_deriv(f);
    2455          91 :   mask = quadratic_prec_mask(e);
    2456          91 :   av2 = avma;
    2457         539 :   for (;mask>1;)
    2458             :   {
    2459             :     GEN u, fa, fr;
    2460         357 :     long n2 = n, rt;
    2461         357 :     n<<=1; if (mask & 1) n--;
    2462         357 :     mask >>= 1;
    2463         357 :     fr = RgXn_red_shallow(f, n);
    2464         357 :     rt = brent_kung_optpow(degpol(fr), 4, 3);
    2465         357 :     an = RgXn_powers(a, rt, n);
    2466         357 :     if (n>1)
    2467             :     {
    2468         357 :       long n4 = (n2+1)>>1;
    2469         357 :       GEN dfr = RgXn_red_shallow(df, n2);
    2470         357 :       dfr = RgX_RgXnV_eval(dfr, RgXnV_red_shallow(an, n2), n2);
    2471         357 :       u = RgX_shift(RgX_Rg_sub(RgXn_mul(W, dfr, n2), gen_1), -n4);
    2472         357 :       W = RgX_sub(W, RgX_shift(RgXn_mul(u, W, n2-n4), n4));
    2473             :     }
    2474         357 :     fa = RgX_sub(RgX_RgXnV_eval(fr, an, n), pol_x(v));
    2475         357 :     fa = RgX_shift(fa, -n2);
    2476         357 :     a = RgX_sub(a, RgX_shift(RgXn_mul(W, fa, n-n2), n2));
    2477         357 :     if (gc_needed(av2,2))
    2478             :     {
    2479           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_reverse, e = %ld", n);
    2480           0 :       gerepileall(av2, 2, &a, &W);
    2481             :     }
    2482             :   }
    2483          91 :   return gerepileupto(av, a);
    2484             : }
    2485             : 
    2486             : GEN
    2487           0 : RgXn_sqrt(GEN h, long e)
    2488             : {
    2489           0 :   pari_sp av = avma, av2;
    2490           0 :   long v = varn(h), n = 1;
    2491           0 :   GEN f = scalarpol(gen_1, v), df = f;
    2492           0 :   ulong mask = quadratic_prec_mask(e);
    2493           0 :   if (degpol(h)<0 || !gequal1(gel(h,2)))
    2494           0 :     pari_err_SQRTN("RgXn_sqrt",h);
    2495           0 :   av2 = avma;
    2496             :   while(1)
    2497           0 :   {
    2498           0 :     long n2 = n, m;
    2499             :     GEN g;
    2500           0 :     n<<=1; if (mask & 1) n--;
    2501           0 :     mask >>= 1;
    2502           0 :     m = n-n2;
    2503           0 :     g = RgX_sub(RgXn_sqrhigh(f, n2, n), RgX_shift_shallow(RgXn_red_shallow(h, n),-n2));
    2504           0 :     f = RgX_sub(f, RgX_shift_shallow(RgXn_mul(gmul2n(df, -1), g, m), n2));
    2505           0 :     if (mask==1) return gerepileupto(av, f);
    2506           0 :     g = RgXn_mul(df, RgXn_mulhigh(df, f, n2, n), m);
    2507           0 :     df = RgX_sub(df, RgX_shift_shallow(g, n2));
    2508           0 :     if (gc_needed(av2,2))
    2509             :     {
    2510           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_sqrt, e = %ld", n);
    2511           0 :       gerepileall(av2, 2, &f, &df);
    2512             :     }
    2513             :   }
    2514             : }
    2515             : 
    2516             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2517             : GEN
    2518       18692 : RgXQ_powu(GEN x, ulong n, GEN T)
    2519             : {
    2520       18692 :   pari_sp av = avma;
    2521             : 
    2522       18692 :   if (!n) return pol_1(varn(x));
    2523       16969 :   if (n == 1) return RgX_copy(x);
    2524       15148 :   x = gen_powu_i(x, n, (void*)T, &_sqr, &_mul);
    2525       15148 :   return gerepilecopy(av, x);
    2526             : }
    2527             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2528             : GEN
    2529       18815 : RgXQ_pow(GEN x, GEN n, GEN T)
    2530             : {
    2531             :   pari_sp av;
    2532       18815 :   long s = signe(n);
    2533             : 
    2534       18815 :   if (!s) return pol_1(varn(x));
    2535       18815 :   if (is_pm1(n) == 1)
    2536           0 :     return (s < 0)? RgXQ_inv(x, T): RgX_copy(x);
    2537       18815 :   av = avma;
    2538       18815 :   if (s < 0) x = RgXQ_inv(x, T);
    2539       18815 :   x = gen_pow_i(x, n, (void*)T, &_sqr, &_mul);
    2540       18815 :   return gerepilecopy(av, x);
    2541             : }
    2542             : static GEN
    2543      159036 : _ZXQsqr(void *data, GEN x) { return ZXQ_sqr(x, (GEN)data); }
    2544             : static GEN
    2545       50465 : _ZXQmul(void *data, GEN x, GEN y) { return ZXQ_mul(x,y, (GEN)data); }
    2546             : /* x,T in Z[X], n in N, compute lift(x^n mod T)) */
    2547             : GEN
    2548      163657 : ZXQ_powu(GEN x, ulong n, GEN T)
    2549             : {
    2550      163657 :   pari_sp av = avma;
    2551             : 
    2552      163657 :   if (!n) return pol_1(varn(x));
    2553      163657 :   if (n == 1) return ZX_copy(x);
    2554      107402 :   x = gen_powu_i(x, n, (void*)T, &_ZXQsqr, &_ZXQmul);
    2555      107408 :   return gerepilecopy(av, x);
    2556             : }
    2557             : 
    2558             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2559             : GEN
    2560        4781 : RgXQ_powers(GEN x, long l, GEN T)
    2561             : {
    2562        4781 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2563        4781 :   return gen_powers(x, l, use_sqr, (void *)T,_sqr,_mul,_one);
    2564             : }
    2565             : 
    2566             : /* a in K = Q[X]/(T), returns [a^0, ..., a^n] */
    2567             : GEN
    2568        3353 : QXQ_powers(GEN a, long n, GEN T)
    2569             : {
    2570        3353 :   GEN den, v = RgXQ_powers(Q_remove_denom(a, &den), n, T);
    2571             :   /* den*a integral; v[i+1] = (den*a)^i in K */
    2572        3353 :   if (den)
    2573             :   { /* restore denominators */
    2574        2044 :     GEN d = den;
    2575             :     long i;
    2576        2044 :     gel(v,2) = a;
    2577        6419 :     for (i=3; i<=n+1; i++) {
    2578        4375 :       d = mulii(d,den);
    2579        4375 :       gel(v,i) = RgX_Rg_div(gel(v,i), d);
    2580             :     }
    2581             :   }
    2582        3353 :   return v;
    2583             : }
    2584             : 
    2585             : static GEN
    2586        2009 : do_QXQ_eval(GEN v, long imin, GEN a, GEN T)
    2587             : {
    2588        2009 :   long l, i, m = 0;
    2589             :   GEN dz, z;
    2590        2009 :   GEN V = cgetg_copy(v, &l);
    2591        6482 :   for (i = imin; i < l; i++)
    2592             :   {
    2593        4473 :     GEN c = gel(v, i);
    2594        4473 :     if (typ(c) == t_POL) m = maxss(m, degpol(c));
    2595             :   }
    2596        2009 :   z = Q_remove_denom(QXQ_powers(a, m, T), &dz);
    2597        2009 :   for (i = 1; i < imin; i++) V[i] = v[i];
    2598        6482 :   for (i = imin; i < l; i++)
    2599             :   {
    2600        4473 :     GEN c = gel(v,i);
    2601        4473 :     if (typ(c) == t_POL) c = QX_ZXQV_eval(c, z, dz);
    2602        4473 :     gel(V,i) = c;
    2603             :   }
    2604        2009 :   return V;
    2605             : }
    2606             : /* [ s(a mod T) | s <- lift(v) ], a,T are QX, v a QXV */
    2607             : GEN
    2608        1946 : QXV_QXQ_eval(GEN v, GEN a, GEN T)
    2609        1946 : { return do_QXQ_eval(v, 1, a, T); }
    2610             : GEN
    2611          63 : QXX_QXQ_eval(GEN v, GEN a, GEN T)
    2612          63 : { return normalizepol(do_QXQ_eval(v, 2, a, T)); }
    2613             : 
    2614             : GEN
    2615         539 : RgXQ_matrix_pow(GEN y, long n, long m, GEN P)
    2616             : {
    2617         539 :   return RgXV_to_RgM(RgXQ_powers(y,m-1,P),n);
    2618             : }
    2619             : 
    2620             : GEN
    2621       42621 : RgXQ_norm(GEN x, GEN T)
    2622             : {
    2623             :   pari_sp av;
    2624       42621 :   long dx = degpol(x);
    2625             :   GEN L, y;
    2626             : 
    2627       42621 :   av = avma; y = resultant(T, x);
    2628       42621 :   L = leading_coeff(T);
    2629       42621 :   if (gequal1(L) || !signe(x)) return y;
    2630           0 :   return gerepileupto(av, gdiv(y, gpowgs(L, dx)));
    2631             : }
    2632             : 
    2633             : GEN
    2634      114925 : RgX_blocks(GEN P, long n, long m)
    2635             : {
    2636      114925 :   GEN z = cgetg(m+1,t_VEC);
    2637      114925 :   long i,j, k=2, l = lg(P);
    2638      537947 :   for(i=1; i<=m; i++)
    2639             :   {
    2640      423022 :     GEN zi = cgetg(n+2,t_POL);
    2641      423022 :     zi[1] = P[1];
    2642      423022 :     gel(z,i) = zi;
    2643     2724865 :     for(j=2; j<n+2; j++)
    2644     2301843 :       gel(zi, j) = k==l ? gen_0 : gel(P,k++);
    2645      423022 :     zi = RgX_renormalize_lg(zi, n+2);
    2646             :   }
    2647      114925 :   return z;
    2648             : }
    2649             : 
    2650             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
    2651             : void
    2652       85322 : RgX_even_odd(GEN p, GEN *pe, GEN *po)
    2653             : {
    2654       85322 :   long n = degpol(p), v = varn(p), n0, n1, i;
    2655             :   GEN p0, p1;
    2656             : 
    2657       85325 :   if (n <= 0) { *pe = RgX_copy(p); *po = zeropol(v); return; }
    2658             : 
    2659       85038 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
    2660       85038 :   p0 = cgetg(n0+2, t_POL); p0[1] = evalvarn(v)|evalsigne(1);
    2661       85040 :   p1 = cgetg(n1+2, t_POL); p1[1] = evalvarn(v)|evalsigne(1);
    2662     5574549 :   for (i=0; i<n1; i++)
    2663             :   {
    2664     5489509 :     p0[2+i] = p[2+(i<<1)];
    2665     5489509 :     p1[2+i] = p[3+(i<<1)];
    2666             :   }
    2667       85040 :   if (n1 != n0)
    2668       45019 :     p0[2+i] = p[2+(i<<1)];
    2669       85040 :   *pe = normalizepol(p0);
    2670       85039 :   *po = normalizepol(p1);
    2671             : }
    2672             : 
    2673             : /* write p(X) = a_0(X^k) + Xa_1(X^k) + ... + X^(k-1)a_{k-1}(X^k), shallow function */
    2674             : GEN
    2675       40922 : RgX_splitting(GEN p, long k)
    2676             : {
    2677       40922 :   long n = degpol(p), v = varn(p), m, i, j, l;
    2678             :   GEN r;
    2679             : 
    2680       40922 :   m = n/k;
    2681       40922 :   r = cgetg(k+1,t_VEC);
    2682      225666 :   for(i=1; i<=k; i++)
    2683             :   {
    2684      184744 :     gel(r,i) = cgetg(m+3, t_POL);
    2685      184744 :     mael(r,i,1) = evalvarn(v)|evalsigne(1);
    2686             :   }
    2687      554799 :   for (j=1, i=0, l=2; i<=n; i++)
    2688             :   {
    2689      513877 :     gmael(r,j,l) = gel(p,2+i);
    2690      513877 :     if (j==k) { j=1; l++; } else j++;
    2691             :   }
    2692      225666 :   for(i=1; i<=k; i++)
    2693      184744 :     gel(r,i) = normalizepol_lg(gel(r,i),i<j?l+1:l);
    2694       40922 :   return r;
    2695             : }
    2696             : 
    2697             : /*******************************************************************/
    2698             : /*                                                                 */
    2699             : /*                        Kronecker form                           */
    2700             : /*                                                                 */
    2701             : /*******************************************************************/
    2702             : 
    2703             : /* z in R[Y] representing an elt in R[X,Y] mod T(Y) in Kronecker form,
    2704             :  * i.e subst(lift(z), x, y^(2deg(z)-1)). Recover the "real" z, with
    2705             :  * normalized coefficients */
    2706             : GEN
    2707      152796 : Kronecker_to_mod(GEN z, GEN T)
    2708             : {
    2709      152796 :   long i,j,lx,l = lg(z), N = (degpol(T)<<1) + 1;
    2710      152796 :   GEN x, t = cgetg(N,t_POL);
    2711      152796 :   t[1] = T[1];
    2712      152796 :   lx = (l-2) / (N-2); x = cgetg(lx+3,t_POL);
    2713      152796 :   x[1] = z[1];
    2714      152796 :   T = RgX_copy(T);
    2715     1485714 :   for (i=2; i<lx+2; i++, z+= N-2)
    2716             :   {
    2717     1332918 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
    2718     1332918 :     gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2719             :   }
    2720      152796 :   N = (l-2) % (N-2) + 2;
    2721      152796 :   for (j=2; j<N; j++) t[j] = z[j];
    2722      152796 :   gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2723      152796 :   return normalizepol_lg(x, i+1);
    2724             : }
    2725             : 
    2726             : /*******************************************************************/
    2727             : /*                                                                 */
    2728             : /*                        Domain detection                         */
    2729             : /*                                                                 */
    2730             : /*******************************************************************/
    2731             : 
    2732             : static GEN
    2733      151257 : zero_FpX_mod(GEN p, long v)
    2734             : {
    2735      151257 :   GEN r = cgetg(3,t_POL);
    2736      151257 :   r[1] = evalvarn(v);
    2737      151257 :   gel(r,2) = mkintmod(gen_0, icopy(p));
    2738      151257 :   return r;
    2739             : }
    2740             : 
    2741             : static GEN
    2742      498854 : RgX_mul_FpX(GEN x, GEN y, GEN p)
    2743             : {
    2744      498854 :   pari_sp av = avma;
    2745             :   GEN r;
    2746      498854 :   if (lgefint(p) == 3)
    2747             :   {
    2748      467784 :     ulong pp = uel(p, 2);
    2749      467784 :     r = Flx_to_ZX_inplace(Flx_mul(RgX_to_Flx(x, pp),
    2750             :                                   RgX_to_Flx(y, pp), pp));
    2751             :   }
    2752             :   else
    2753       31070 :     r = FpX_mul(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2754      498854 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2755      385922 :   return gerepileupto(av, FpX_to_mod(r, p));
    2756             : }
    2757             : 
    2758             : static GEN
    2759         532 : zero_FpXQX_mod(GEN pol, GEN p, long v)
    2760             : {
    2761         532 :   GEN r = cgetg(3,t_POL);
    2762         532 :   r[1] = evalvarn(v);
    2763         532 :   gel(r,2) = mkpolmod(mkintmod(gen_0, icopy(p)), gcopy(pol));
    2764         532 :   return r;
    2765             : }
    2766             : 
    2767             : static GEN
    2768        2037 : RgX_mul_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2769             : {
    2770        2037 :   pari_sp av = avma;
    2771             :   long dT;
    2772             :   GEN kx, ky, r;
    2773        2037 :   GEN T = RgX_to_FpX(pol, p);
    2774        2037 :   if (signe(T)==0) pari_err_OP("*", x, y);
    2775        2030 :   dT = degpol(T);
    2776        2030 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2777        2030 :   ky = ZXX_to_Kronecker(RgX_to_FpXQX(y, T, p), dT);
    2778        2030 :   r = FpX_mul(kx, ky, p);
    2779        2030 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2780        1505 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2781             : }
    2782             : 
    2783             : static GEN
    2784      369687 : RgX_liftred(GEN x, GEN T)
    2785      369687 : { return RgXQX_red(liftpol_shallow(x), T); }
    2786             : 
    2787             : static GEN
    2788      138773 : RgX_mul_QXQX(GEN x, GEN y, GEN T)
    2789             : {
    2790      138773 :   pari_sp av = avma;
    2791      138773 :   long dT = degpol(T);
    2792      138773 :   GEN r = QX_mul(ZXX_to_Kronecker(RgX_liftred(x, T), dT),
    2793             :                  ZXX_to_Kronecker(RgX_liftred(y, T), dT));
    2794      138773 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2795             : }
    2796             : 
    2797             : static GEN
    2798        1260 : RgX_sqr_FpX(GEN x, GEN p)
    2799             : {
    2800        1260 :   pari_sp av = avma;
    2801             :   GEN r;
    2802        1260 :   if (lgefint(p) == 3)
    2803             :   {
    2804        1056 :     ulong pp = uel(p, 2);
    2805        1056 :     r = Flx_to_ZX_inplace(Flx_sqr(RgX_to_Flx(x, pp), pp));
    2806             :   }
    2807             :   else
    2808         204 :     r = FpX_sqr(RgX_to_FpX(x, p), p);
    2809        1260 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2810        1162 :   return gerepileupto(av, FpX_to_mod(r, p));
    2811             : }
    2812             : 
    2813             : static GEN
    2814         196 : RgX_sqr_FpXQX(GEN x, GEN pol, GEN p)
    2815             : {
    2816         196 :   pari_sp av = avma;
    2817             :   long dT;
    2818         196 :   GEN kx, r, T = RgX_to_FpX(pol, p);
    2819         196 :   if (signe(T)==0) pari_err_OP("*",x,x);
    2820         189 :   dT = degpol(T);
    2821         189 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2822         189 :   r = FpX_sqr(kx, p);
    2823         189 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2824         189 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2825             : }
    2826             : 
    2827             : static GEN
    2828       12329 : RgX_sqr_QXQX(GEN x, GEN T)
    2829             : {
    2830       12329 :   pari_sp av = avma;
    2831       12329 :   long dT = degpol(T);
    2832       12329 :   GEN r = QX_sqr(ZXX_to_Kronecker(RgX_liftred(x, T), dT));
    2833       12329 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2834             : }
    2835             : 
    2836             : static GEN
    2837       69951 : RgX_rem_FpX(GEN x, GEN y, GEN p)
    2838             : {
    2839       69951 :   pari_sp av = avma;
    2840             :   GEN r;
    2841       69951 :   if (lgefint(p) == 3)
    2842             :   {
    2843       50690 :     ulong pp = uel(p, 2);
    2844       50690 :     r = Flx_to_ZX_inplace(Flx_rem(RgX_to_Flx(x, pp),
    2845             :                                   RgX_to_Flx(y, pp), pp));
    2846             :   }
    2847             :   else
    2848       19261 :     r = FpX_rem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2849       69951 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2850       31724 :   return gerepileupto(av, FpX_to_mod(r, p));
    2851             : }
    2852             : 
    2853             : static GEN
    2854       39906 : RgX_rem_QXQX(GEN x, GEN y, GEN T)
    2855             : {
    2856       39906 :   pari_sp av = avma;
    2857             :   GEN r;
    2858       39906 :   r = RgXQX_rem(RgX_liftred(x, T), RgX_liftred(y, T), T);
    2859       39906 :   return gerepilecopy(av, QXQX_to_mod_shallow(r, T));
    2860             : }
    2861             : static GEN
    2862          63 : RgX_rem_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2863             : {
    2864          63 :   pari_sp av = avma;
    2865             :   GEN r;
    2866          63 :   GEN T = RgX_to_FpX(pol, p);
    2867          63 :   if (signe(T) == 0) pari_err_OP("%", x, y);
    2868          56 :   if (lgefint(p) == 3)
    2869             :   {
    2870          48 :     ulong pp = uel(p, 2);
    2871          48 :     GEN Tp = ZX_to_Flx(T, pp);
    2872          48 :     r = FlxX_to_ZXX(FlxqX_rem(RgX_to_FlxqX(x, Tp, pp),
    2873             :                               RgX_to_FlxqX(y, Tp, pp), Tp, pp));
    2874             :   }
    2875             :   else
    2876           8 :     r = FpXQX_rem(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), T, p);
    2877          56 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2878          49 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2879             : }
    2880             : 
    2881             : #define code(t1,t2) ((t1 << 6) | t2)
    2882             : static GEN
    2883    48541445 : RgX_mul_fast(GEN x, GEN y)
    2884             : {
    2885             :   GEN p, pol;
    2886             :   long pa;
    2887    48541445 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2888    48541443 :   switch(t)
    2889             :   {
    2890    36350460 :     case t_INT:    return ZX_mul(x,y);
    2891     1571928 :     case t_FRAC:   return QX_mul(x,y);
    2892      104161 :     case t_FFELT:  return FFX_mul(x, y, pol);
    2893      498854 :     case t_INTMOD: return RgX_mul_FpX(x, y, p);
    2894             :     case code(t_POLMOD, t_INT):
    2895             :     case code(t_POLMOD, t_FRAC):
    2896      138773 :                    return RgX_mul_QXQX(x, y, pol);
    2897             :     case code(t_POLMOD, t_INTMOD):
    2898        2037 :                    return RgX_mul_FpXQX(x, y, pol, p);
    2899     9875230 :     default:       return NULL;
    2900             :   }
    2901             : }
    2902             : static GEN
    2903      143760 : RgX_sqr_fast(GEN x)
    2904             : {
    2905             :   GEN p, pol;
    2906             :   long pa;
    2907      143760 :   long t = RgX_type(x,&p,&pol,&pa);
    2908      143760 :   switch(t)
    2909             :   {
    2910      109759 :     case t_INT:    return ZX_sqr(x);
    2911       11375 :     case t_FRAC:   return QX_sqr(x);
    2912        2296 :     case t_FFELT:  return FFX_sqr(x, pol);
    2913        1260 :     case t_INTMOD: return RgX_sqr_FpX(x, p);
    2914             :     case code(t_POLMOD, t_INT):
    2915             :     case code(t_POLMOD, t_FRAC):
    2916       12329 :                    return RgX_sqr_QXQX(x, pol);
    2917             :     case code(t_POLMOD, t_INTMOD):
    2918         196 :                    return RgX_sqr_FpXQX(x, pol, p);
    2919        6545 :     default:       return NULL;
    2920             :   }
    2921             : }
    2922             : 
    2923             : static GEN
    2924     6307618 : RgX_rem_fast(GEN x, GEN y)
    2925             : {
    2926             :   GEN p, pol;
    2927             :   long pa;
    2928     6307618 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2929     6307618 :   switch(t)
    2930             :   {
    2931     3933732 :     case t_INT:    return ZX_is_monic(y) ? ZX_rem(x,y): NULL;
    2932     1231768 :     case t_FRAC:   return RgX_is_ZX(y) && ZX_is_monic(y) ? QX_ZX_rem(x,y): NULL;
    2933          84 :     case t_FFELT:  return FFX_rem(x, y, pol);
    2934       69951 :     case t_INTMOD: return RgX_rem_FpX(x, y, p);
    2935             :     case code(t_POLMOD, t_INT):
    2936             :     case code(t_POLMOD, t_FRAC):
    2937       39906 :                    return RgX_rem_QXQX(x, y, pol);
    2938             :     case code(t_POLMOD, t_INTMOD):
    2939          63 :                    return RgX_rem_FpXQX(x, y, pol, p);
    2940     1032114 :     default:       return NULL;
    2941             :   }
    2942             : }
    2943             : 
    2944             : #undef code
    2945             : 
    2946             : GEN
    2947    46068989 : RgX_mul(GEN x, GEN y)
    2948             : {
    2949    46068989 :   GEN z = RgX_mul_fast(x,y);
    2950    46068979 :   if (!z) z = RgX_mul_i(x,y);
    2951    46068979 :   return z;
    2952             : }
    2953             : 
    2954             : GEN
    2955      127093 : RgX_sqr(GEN x)
    2956             : {
    2957      127093 :   GEN z = RgX_sqr_fast(x);
    2958      127086 :   if (!z) z = RgX_sqr_i(x);
    2959      127086 :   return z;
    2960             : }
    2961             : 
    2962             : GEN
    2963     6307618 : RgX_rem(GEN x, GEN y)
    2964             : {
    2965     6307618 :   GEN z = RgX_rem_fast(x, y);
    2966     6307611 :   if (!z) z = RgX_divrem_i(x, y, ONLY_REM);
    2967     6307611 :   return z;
    2968             : }
    2969             : 
    2970             : GEN
    2971     2427643 : RgXn_mul(GEN f, GEN g, long n)
    2972             : {
    2973     2427643 :   pari_sp av = avma;
    2974     2427643 :   GEN h = RgX_mul_fast(f,g);
    2975     2427643 :   if (!h) return RgXn_mul2(f,g,n);
    2976      173545 :   if (degpol(h) < n) return h;
    2977       48986 :   return gerepilecopy(av, RgXn_red_shallow(h, n));
    2978             : }
    2979             : 
    2980             : GEN
    2981       16667 : RgXn_sqr(GEN f, long n)
    2982             : {
    2983       16667 :   pari_sp av = avma;
    2984       16667 :   GEN g = RgX_sqr_fast(f);
    2985       16667 :   if (!g) return RgXn_sqr2(f,n);
    2986       16429 :   if (degpol(g) < n) return g;
    2987        1085 :   return gerepilecopy(av, RgXn_red_shallow(g, n));
    2988             : }
    2989             : 
    2990             : /* (f*g) \/ x^n */
    2991             : GEN
    2992       44813 : RgX_mulhigh_i(GEN f, GEN g, long n)
    2993             : {
    2994       44813 :   GEN h = RgX_mul_fast(f,g);
    2995       44813 :   return h? RgX_shift_shallow(h, -n): RgX_mulhigh_i2(f,g,n);
    2996             : }
    2997             : 
    2998             : /* (f*g) \/ x^n */
    2999             : GEN
    3000           0 : RgX_sqrhigh_i(GEN f, long n)
    3001             : {
    3002           0 :   GEN h = RgX_sqr_fast(f);
    3003           0 :   return h? RgX_shift_shallow(h, -n): RgX_sqrhigh_i2(f,n);
    3004             : }

Generated by: LCOV version 1.13