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.0 lcov report (development 23036-b751c0af5) Lines: 1460 1627 89.7 %
Date: 2018-09-26 05:46:06 Functions: 175 190 92.1 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8             : ANY WARRANTY WHATSOEVER.
       9             : 
      10             : Check the License for details. You should have received a copy of it, along
      11             : with the package; see the file 'COPYING'. If not, write to the Free Software
      12             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13             : 
      14             : #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     3175670 : brent_kung_optpow(long d, long n, long m)
      29             : {
      30             :   long p, r;
      31     3175670 :   long pold=1, rold=n*(d-1);
      32    18481975 :   for(p=2; p<=d; p++)
      33             :   {
      34    15306305 :     r = m*(p-1) + n*((d-1)/p);
      35    15306305 :     if (r<rold) { pold=p; rold=r; }
      36             :   }
      37     3175670 :   return pold;
      38             : }
      39             : 
      40             : static GEN
      41     7064846 : 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     7064846 :   pari_sp av = avma;
      45             :   long i;
      46     7064846 :   GEN z = cmul(E,P,a,ff->one(E));
      47     7065446 :   if (!z) z = gen_0;
      48    50256052 :   for (i=1; i<=n; i++)
      49             :   {
      50    43191247 :     GEN t = cmul(E,P,a+i,gel(V,i+1));
      51    43194402 :     if (t) {
      52    41402700 :       z = ff->add(E, z, t);
      53    41398890 :       if (gc_needed(av,2)) z = gerepileupto(av, z);
      54             :     }
      55             :   }
      56     7064805 :   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     5866771 : 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     5866771 :   pari_sp av = avma;
      69     5866771 :   long l = lg(V)-1;
      70             :   GEN z, u;
      71             : 
      72     5866771 :   if (d < 0) return ff->zero(E);
      73     5557780 :   if (d < l) return gerepileupto(av, gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul));
      74      695188 :   if (l<2) pari_err_DOMAIN("gen_RgX_bkeval_powers", "#powers", "<",gen_2,V);
      75      695188 :   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      695188 :   d -= l;
      81      695188 :   z = gen_RgXQ_eval_powers(P,V,d+1,l-1,E,ff,cmul);
      82     2202324 :   while (d >= l-1)
      83             :   {
      84      811953 :     d -= l-1;
      85      811953 :     u = gen_RgXQ_eval_powers(P,V,d+1,l-2,E,ff,cmul);
      86      811895 :     z = ff->add(E,u, ff->mul(E,z,gel(V,l)));
      87      811939 :     if (gc_needed(av,2))
      88          61 :       z = gerepileupto(av, z);
      89             :   }
      90      695188 :   u = gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul);
      91      695190 :   z = ff->add(E,u, ff->mul(E,z,gel(V,d+2)));
      92      695188 :   return gerepileupto(av, ff->red(E,z));
      93             : }
      94             : 
      95             : GEN
      96      450562 : 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      450562 :   pari_sp av = avma;
     100             :   GEN z, V;
     101             :   long rtd;
     102      450562 :   if (d < 0) return ff->zero(E);
     103      450212 :   rtd = (long) sqrt((double)d);
     104      450212 :   V = gen_powers(x,rtd,use_sqr,E,ff->sqr,ff->mul,ff->one);
     105      450206 :   z = gen_bkeval_powers(Q, d, V, E, ff, cmul);
     106      450210 :   return gerepileupto(av, z);
     107             : }
     108             : 
     109             : static GEN
     110      903605 : _gen_nored(void *E, GEN x) { (void)E; return x; }
     111             : static GEN
     112    21775627 : _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      762435 : _gen_mul(void *E, GEN x, GEN y) { (void)E; return gmul(x, y); }
     117             : static GEN
     118      250294 : _gen_sqr(void *E, GEN x) { (void)E; return gsqr(x); }
     119             : static GEN
     120      921455 : _gen_one(void *E) { (void)E; return gen_1; }
     121             : static GEN
     122       14735 : _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      471912 : _gen_cmul(void *E, GEN P, long a, GEN x)
     129      471912 : {(void)E; return gmul(gel(P,a+2), x);}
     130             : 
     131             : GEN
     132      143465 : RgX_RgV_eval(GEN Q, GEN x)
     133             : {
     134      143465 :   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        2884 : RgXV_RgV_eval(GEN Q, GEN x)
     145             : {
     146        2884 :   long i, l = lg(Q), vQ = gvar(Q);
     147        2884 :   GEN v = cgetg(l, t_VEC);
     148      242928 :   for (i = 1; i < l; i++)
     149             :   {
     150      240044 :     GEN Qi = gel(Q, i);
     151      240044 :     gel(v, i) = typ(Qi)==t_POL && varn(Qi)==vQ? RgX_RgV_eval(Qi, x): gcopy(Qi);
     152             :   }
     153        2884 :   return v;
     154             : }
     155             : 
     156             : const struct bb_algebra *
     157      115175 : get_Rg_algebra(void)
     158             : {
     159      115175 :   return &Rg_algebra;
     160             : }
     161             : 
     162             : static struct bb_ring Rg_ring = {  _gen_add, _gen_mul, _gen_sqr };
     163             : 
     164             : static GEN
     165       17668 : _RgX_divrem(void *E, GEN x, GEN y, GEN *r)
     166             : {
     167             :   (void) E;
     168       17668 :   return RgX_divrem(x, y, r);
     169             : }
     170             : 
     171             : GEN
     172        4221 : RgX_digits(GEN x, GEN T)
     173             : {
     174        4221 :   pari_sp av = avma;
     175        4221 :   long d = degpol(T), n = (lgpol(x)+d-1)/d;
     176        4221 :   GEN z = gen_digits(x,T,n,NULL, &Rg_ring, _RgX_divrem);
     177        4221 :   return gerepileupto(av, z);
     178             : }
     179             : 
     180             : /*******************************************************************/
     181             : /*                                                                 */
     182             : /*                         RgX                                     */
     183             : /*                                                                 */
     184             : /*******************************************************************/
     185             : 
     186             : long
     187     4397510 : RgX_equal(GEN x, GEN y)
     188             : {
     189     4397510 :   long i = lg(x);
     190             : 
     191     4397510 :   if (i != lg(y)) return 0;
     192    22628112 :   for (i--; i > 1; i--)
     193    18357247 :     if (!gequal(gel(x,i),gel(y,i))) return 0;
     194     4270865 :   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     2384648 : Rg_get_1(GEN x)
     201             : {
     202             :   GEN p, T;
     203     2384648 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     204     2384648 :   if (RgX_type_is_composite(tx))
     205      349832 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     206     2384648 :   switch(tx)
     207             :   {
     208         161 :     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         637 :     case t_FFELT: return FF_1(T);
     211     2383843 :     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     7032312 : Rg_get_0(GEN x)
     218             : {
     219             :   GEN p, T;
     220     7032312 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     221     7032312 :   if (RgX_type_is_composite(tx))
     222       75166 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     223     7032312 :   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     7031759 :     default: return gen_0;
     229             :   }
     230             : }
     231             : 
     232             : GEN
     233        3311 : QX_ZXQV_eval(GEN P, GEN V, GEN dV)
     234             : {
     235        3311 :   long i, n = degpol(P);
     236             :   GEN z, dz, dP;
     237        3311 :   if (n < 0) return gen_0;
     238        3311 :   P = Q_remove_denom(P, &dP);
     239        3311 :   z = gel(P,2); if (n == 0) return icopy(z);
     240        1778 :   if (dV) z = mulii(dV, z); /* V[1] = dV */
     241        1778 :   z = ZX_Z_add_shallow(ZX_Z_mul(gel(V,2),gel(P,3)), z);
     242        1778 :   for (i=2; i<=n; i++) z = ZX_add(ZX_Z_mul(gel(V,i+1),gel(P,2+i)), z);
     243        1778 :   dz = mul_denom(dP, dV);
     244        1778 :   return dz? RgX_Rg_div(z, dz): z;
     245             : }
     246             : 
     247             : /* Return P(h * x), not memory clean */
     248             : GEN
     249       10486 : RgX_unscale(GEN P, GEN h)
     250             : {
     251       10486 :   long i, l = lg(P);
     252       10486 :   GEN hi = gen_1, Q = cgetg(l, t_POL);
     253       10486 :   Q[1] = P[1];
     254       10486 :   if (l == 2) return Q;
     255       10458 :   gel(Q,2) = gcopy(gel(P,2));
     256       35448 :   for (i=3; i<l; i++)
     257             :   {
     258       24990 :     hi = gmul(hi,h);
     259       24990 :     gel(Q,i) = gmul(gel(P,i), hi);
     260             :   }
     261       10458 :   return Q;
     262             : }
     263             : /* P a ZX, Return P(h * x), not memory clean; optimize for h = -1 */
     264             : GEN
     265      359688 : ZX_z_unscale(GEN P, long h)
     266             : {
     267      359688 :   long i, l = lg(P);
     268      359688 :   GEN Q = cgetg(l, t_POL);
     269      359688 :   Q[1] = P[1];
     270      359688 :   if (l == 2) return Q;
     271      359688 :   gel(Q,2) = gel(P,2);
     272      359688 :   if (l == 3) return Q;
     273      359688 :   if (h == -1)
     274      309449 :     for (i = 3; i < l; i++)
     275             :     {
     276      301391 :       gel(Q,i) = negi(gel(P,i));
     277      301391 :       if (++i == l) break;
     278      298011 :       gel(Q,i) = gel(P,i);
     279             :     }
     280             :   else
     281             :   {
     282             :     GEN hi;
     283      348250 :     gel(Q,3) = mulis(gel(P,3), h);
     284      348250 :     hi = sqrs(h);
     285     1704808 :     for (i = 4; i < l; i++)
     286             :     {
     287     1356558 :       gel(Q,i) = mulii(gel(P,i), hi);
     288     1356558 :       if (i != l-1) hi = mulis(hi,h);
     289             :     }
     290             :   }
     291      359688 :   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        7343 : ZX_unscale(GEN P, GEN h)
     296             : {
     297             :   long i, l;
     298             :   GEN Q, hi;
     299        7343 :   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       56806 : ZX_unscale2n(GEN P, long n)
     317             : {
     318       56806 :   long i, ni = n, l = lg(P);
     319       56806 :   GEN Q = cgetg(l, t_POL);
     320       56806 :   Q[1] = P[1];
     321       56806 :   if (l == 2) return Q;
     322       56806 :   gel(Q,2) = gel(P,2);
     323       56806 :   if (l == 3) return Q;
     324       56806 :   gel(Q,3) = shifti(gel(P,3), ni);
     325      226349 :   for (i=4; i<l; i++)
     326             :   {
     327      169543 :     ni += n;
     328      169543 :     gel(Q,i) = shifti(gel(P,i), ni);
     329             :   }
     330       56806 :   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         385 : RgXV_unscale(GEN v, GEN h)
     356             : {
     357             :   long i, l;
     358             :   GEN w;
     359         385 :   if (!h || isint1(h)) return v;
     360         308 :   w = cgetg_copy(v, &l);
     361         308 :   for (i=1; i<l; i++) gel(w,i) = RgX_unscale(gel(v,i), h);
     362         308 :   return w;
     363             : }
     364             : 
     365             : /* Return h^degpol(P) P(x / h), not memory clean */
     366             : GEN
     367        2002 : RgX_rescale(GEN P, GEN h)
     368             : {
     369        2002 :   long i, l = lg(P);
     370        2002 :   GEN Q = cgetg(l,t_POL), hi = h;
     371        2002 :   Q[l-1] = P[l-1];
     372        9142 :   for (i=l-2; i>=2; i--)
     373             :   {
     374        9142 :     gel(Q,i) = gmul(gel(P,i), hi);
     375        9142 :     if (i == 2) break;
     376        7140 :     hi = gmul(hi,h);
     377             :   }
     378        2002 :   Q[1] = P[1]; return Q;
     379             : }
     380             : 
     381             : /* A(X^d) --> A(X) */
     382             : GEN
     383      126518 : RgX_deflate(GEN x0, long d)
     384             : {
     385             :   GEN z, y, x;
     386      126518 :   long i,id, dy, dx = degpol(x0);
     387      126518 :   if (d == 1 || dx <= 0) return leafcopy(x0);
     388       70714 :   dy = dx/d;
     389       70714 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     390       70714 :   z = y + 2;
     391       70714 :   x = x0+ 2;
     392       70714 :   for (i=id=0; i<=dy; i++,id+=d) gel(z,i) = gel(x,id);
     393       70714 :   return y;
     394             : }
     395             : 
     396             : /* return x0(X^d) */
     397             : GEN
     398      458182 : RgX_inflate(GEN x0, long d)
     399             : {
     400      458182 :   long i, id, dy, dx = degpol(x0);
     401      458182 :   GEN x = x0 + 2, z, y;
     402      458182 :   if (dx <= 0) return leafcopy(x0);
     403      451000 :   dy = dx*d;
     404      451000 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     405      451000 :   z = y + 2;
     406      451000 :   for (i=0; i<=dy; i++) gel(z,i) = gen_0;
     407      451000 :   for (i=id=0; i<=dx; i++,id+=d) gel(z,id) = gel(x,i);
     408      451000 :   return y;
     409             : }
     410             : 
     411             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     412             : GEN
     413     1363813 : RgX_translate(GEN P, GEN c)
     414             : {
     415     1363813 :   pari_sp av = avma;
     416             :   GEN Q, *R;
     417             :   long i, k, n;
     418             : 
     419     1363813 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     420     1362158 :   Q = leafcopy(P);
     421     1362158 :   R = (GEN*)(Q+2); n = degpol(P);
     422     1362158 :   if (equali1(c))
     423             :   {
     424        8862 :     for (i=1; i<=n; i++)
     425             :     {
     426        6713 :       for (k=n-i; k<n; k++) R[k] = gadd(R[k], R[k+1]);
     427        6713 :       if (gc_needed(av,2))
     428             :       {
     429           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL(1), i = %ld/%ld", i,n);
     430           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     431             :       }
     432             :     }
     433             :   }
     434     1360009 :   else if (equalim1(c))
     435             :   {
     436      142352 :     for (i=1; i<=n; i++)
     437             :     {
     438      120827 :       for (k=n-i; k<n; k++) R[k] = gsub(R[k], R[k+1]);
     439      120827 :       if (gc_needed(av,2))
     440             :       {
     441           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL(-1), i = %ld/%ld", i,n);
     442           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     443             :       }
     444             :     }
     445             :   }
     446             :   else
     447             :   {
     448     4543840 :     for (i=1; i<=n; i++)
     449             :     {
     450     3205356 :       for (k=n-i; k<n; k++) R[k] = gadd(R[k], gmul(c, R[k+1]));
     451     3205356 :       if (gc_needed(av,2))
     452             :       {
     453           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"TR_POL, i = %ld/%ld", i,n);
     454           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     455             :       }
     456             :     }
     457             :   }
     458     1362158 :   return gerepilecopy(av, Q);
     459             : }
     460             : 
     461             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     462             : GEN
     463      388497 : ZX_translate(GEN P, GEN c)
     464             : {
     465      388497 :   pari_sp av = avma;
     466             :   GEN Q, *R;
     467             :   long i, k, n;
     468             : 
     469      388497 :   if (!signe(P) || !signe(c)) return ZX_copy(P);
     470      388462 :   Q = leafcopy(P);
     471      388462 :   R = (GEN*)(Q+2); n = degpol(P);
     472      388462 :   if (equali1(c))
     473             :   {
     474     2458535 :     for (i=1; i<=n; i++)
     475             :     {
     476     2174291 :       for (k=n-i; k<n; k++) R[k] = addii(R[k], R[k+1]);
     477     2174291 :       if (gc_needed(av,2))
     478             :       {
     479           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(1), i = %ld/%ld", i,n);
     480           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     481             :       }
     482             :     }
     483             :   }
     484      104218 :   else if (equalim1(c))
     485             :   {
     486          70 :     for (i=1; i<=n; i++)
     487             :     {
     488          49 :       for (k=n-i; k<n; k++) R[k] = subii(R[k], R[k+1]);
     489          49 :       if (gc_needed(av,2))
     490             :       {
     491           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate(-1), i = %ld/%ld", i,n);
     492           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     493             :       }
     494             :     }
     495             :   }
     496             :   else
     497             :   {
     498      782615 :     for (i=1; i<=n; i++)
     499             :     {
     500      678418 :       for (k=n-i; k<n; k++) R[k] = addii(R[k], mulii(c, R[k+1]));
     501      678418 :       if (gc_needed(av,2))
     502             :       {
     503           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"ZX_translate, i = %ld/%ld", i,n);
     504           0 :         Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     505             :       }
     506             :     }
     507             :   }
     508      388462 :   return gerepilecopy(av, Q);
     509             : }
     510             : /* return lift( P(X + c) ) using Horner, c in R[y]/(T) */
     511             : GEN
     512        6426 : RgXQX_translate(GEN P, GEN c, GEN T)
     513             : {
     514        6426 :   pari_sp av = avma;
     515             :   GEN Q, *R;
     516             :   long i, k, n;
     517             : 
     518        6426 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     519        6405 :   Q = leafcopy(P);
     520        6405 :   R = (GEN*)(Q+2); n = degpol(P);
     521       38311 :   for (i=1; i<=n; i++)
     522             :   {
     523      164997 :     for (k=n-i; k<n; k++)
     524             :     {
     525      133091 :       pari_sp av2 = avma;
     526      133091 :       R[k] = gerepileupto(av2, RgX_rem(gadd(R[k], gmul(c, R[k+1])), T));
     527             :     }
     528       31906 :     if (gc_needed(av,2))
     529             :     {
     530           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXQX_translate, i = %ld/%ld", i,n);
     531           0 :       Q = gerepilecopy(av, Q); R = (GEN*)Q+2;
     532             :     }
     533             :   }
     534        6405 :   return gerepilecopy(av, Q);
     535             : }
     536             : 
     537             : /********************************************************************/
     538             : /**                                                                **/
     539             : /**                          CONVERSIONS                           **/
     540             : /**                       (not memory clean)                       **/
     541             : /**                                                                **/
     542             : /********************************************************************/
     543             : /* to INT / FRAC / (POLMOD mod T), not memory clean because T not copied,
     544             :  * but everything else is */
     545             : static GEN
     546       24021 : QXQ_to_mod_copy(GEN x, GEN T)
     547             : {
     548             :   long d;
     549       24021 :   switch(typ(x))
     550             :   {
     551        9660 :     case t_INT:  return icopy(x);
     552         805 :     case t_FRAC: return gcopy(x);
     553             :     case t_POL:
     554       13556 :       d = degpol(x);
     555       13556 :       if (d < 0) return gen_0;
     556       13017 :       if (d == 0) return gcopy(gel(x,2));
     557       12667 :       return mkpolmod(RgX_copy(x), T);
     558           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     559             :              return NULL;/* LCOV_EXCL_LINE */
     560             :   }
     561             : }
     562             : /* pure shallow version */
     563             : static GEN
     564      978813 : QXQ_to_mod(GEN x, GEN T)
     565             : {
     566             :   long d;
     567      978813 :   switch(typ(x))
     568             :   {
     569             :     case t_INT:
     570      433489 :     case t_FRAC: return x;
     571             :     case t_POL:
     572      545324 :       d = degpol(x);
     573      545324 :       if (d < 0) return gen_0;
     574      416720 :       if (d == 0) return gel(x,2);
     575      399487 :       return mkpolmod(x, T);
     576           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     577             :              return NULL;/* LCOV_EXCL_LINE */
     578             :   }
     579             : }
     580             : /* T a ZX, z lifted from (Q[Y]/(T(Y)))[X], apply QXQ_to_mod_copy to all coeffs.
     581             :  * Not memory clean because T not copied, but everything else is */
     582             : static GEN
     583        3206 : QXQX_to_mod(GEN z, GEN T)
     584             : {
     585        3206 :   long i,l = lg(z);
     586        3206 :   GEN x = cgetg(l,t_POL);
     587        3206 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod_copy(gel(z,i), T);
     588        3206 :   x[1] = z[1]; return normalizepol_lg(x,l);
     589             : }
     590             : /* pure shallow version */
     591             : GEN
     592      146663 : QXQX_to_mod_shallow(GEN z, GEN T)
     593             : {
     594      146663 :   long i,l = lg(z);
     595      146663 :   GEN x = cgetg(l,t_POL);
     596      146663 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod(gel(z,i), T);
     597      146663 :   x[1] = z[1]; return normalizepol_lg(x,l);
     598             : }
     599             : /* Apply QXQX_to_mod to all entries. Memory-clean ! */
     600             : GEN
     601        1239 : QXQXV_to_mod(GEN V, GEN T)
     602             : {
     603        1239 :   long i, l = lg(V);
     604        1239 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     605        1239 :   for (i=1;i<l; i++) gel(z,i) = QXQX_to_mod(gel(V,i), T);
     606        1239 :   return z;
     607             : }
     608             : /* Apply QXQ_to_mod_copy to all entries. Memory-clean ! */
     609             : GEN
     610        4908 : QXQV_to_mod(GEN V, GEN T)
     611             : {
     612        4908 :   long i, l = lg(V);
     613        4908 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     614        4908 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod_copy(gel(V,i), T);
     615        4908 :   return z;
     616             : }
     617             : 
     618             : /* Apply QXQ_to_mod_copy to all entries. Memory-clean ! */
     619             : GEN
     620       19222 : QXQC_to_mod_shallow(GEN V, GEN T)
     621             : {
     622       19222 :   long i, l = lg(V);
     623       19222 :   GEN z = cgetg(l, t_COL);
     624       19222 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod(gel(V,i), T);
     625       19222 :   return z;
     626             : }
     627             : 
     628             : GEN
     629        6615 : QXQM_to_mod_shallow(GEN V, GEN T)
     630             : {
     631        6615 :   long i, l = lg(V);
     632        6615 :   GEN z = cgetg(l, t_MAT);
     633        6615 :   for (i=1; i<l; i++) gel(z,i) = QXQC_to_mod_shallow(gel(V,i), T);
     634        6615 :   return z;
     635             : }
     636             : 
     637             : GEN
     638      778653 : RgX_renormalize_lg(GEN x, long lx)
     639             : {
     640             :   long i;
     641     2128906 :   for (i = lx-1; i>1; i--)
     642     1991264 :     if (! gequal0(gel(x,i))) break; /* _not_ isexactzero */
     643      778653 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
     644      778653 :   setlg(x, i+1); setsigne(x, i != 1); return x;
     645             : }
     646             : 
     647             : GEN
     648      468203 : RgV_to_RgX(GEN x, long v)
     649             : {
     650      468203 :   long i, k = lg(x);
     651             :   GEN p;
     652             : 
     653      468203 :   while (--k && gequal0(gel(x,k)));
     654      468201 :   if (!k) return pol_0(v);
     655      462482 :   i = k+2; p = cgetg(i,t_POL);
     656      462482 :   p[1] = evalsigne(1) | evalvarn(v);
     657      462482 :   x--; for (k=2; k<i; k++) gel(p,k) = gel(x,k);
     658      462482 :   return p;
     659             : }
     660             : GEN
     661      158603 : RgV_to_RgX_reverse(GEN x, long v)
     662             : {
     663      158603 :   long j, k, l = lg(x);
     664             :   GEN p;
     665             : 
     666      158603 :   for (k = 1; k < l; k++)
     667      158603 :     if (!gequal0(gel(x,k))) break;
     668      158603 :   if (k == l) return pol_0(v);
     669      158603 :   k -= 1;
     670      158603 :   l -= k;
     671      158603 :   x += k;
     672      158603 :   p = cgetg(l+1,t_POL);
     673      158603 :   p[1] = evalsigne(1) | evalvarn(v);
     674      158603 :   for (j=2, k=l; j<=l; j++) gel(p,j) = gel(x,--k);
     675      158603 :   return p;
     676             : }
     677             : 
     678             : /* return the (N-dimensional) vector of coeffs of p */
     679             : GEN
     680     6857821 : RgX_to_RgC(GEN x, long N)
     681             : {
     682             :   long i, l;
     683             :   GEN z;
     684     6857821 :   l = lg(x)-1; x++;
     685     6857821 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
     686     6857821 :   z = cgetg(N+1,t_COL);
     687     6857821 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
     688     6857821 :   for (   ; i<=N; i++) gel(z,i) = gen_0;
     689     6857821 :   return z;
     690             : }
     691             : GEN
     692      689918 : Rg_to_RgC(GEN x, long N)
     693             : {
     694      689918 :   return (typ(x) == t_POL)? RgX_to_RgC(x,N): scalarcol_shallow(x, N);
     695             : }
     696             : 
     697             : /* vector of polynomials (in v) whose coeffs are given by the columns of x */
     698             : GEN
     699       49203 : RgM_to_RgXV(GEN x, long v)
     700       49203 : { pari_APPLY_type(t_VEC, RgV_to_RgX(gel(x,i), v)) }
     701             : 
     702             : /* matrix whose entries are given by the coeffs of the polynomials in
     703             :  * vector v (considered as degree n-1 polynomials) */
     704             : GEN
     705      105777 : RgV_to_RgM(GEN x, long n)
     706      105777 : { pari_APPLY_type(t_MAT, Rg_to_RgC(gel(x,i), n)) }
     707             : 
     708             : GEN
     709        3984 : RgXV_to_RgM(GEN x, long n)
     710        3984 : { pari_APPLY_type(t_MAT, RgX_to_RgC(gel(x,i), n)) }
     711             : 
     712             : /* polynomial (in v) of polynomials (in w) whose coeffs are given by the columns of x */
     713             : GEN
     714       16359 : RgM_to_RgXX(GEN x, long v,long w)
     715             : {
     716       16359 :   long j, lx = lg(x);
     717       16359 :   GEN y = cgetg(lx+1, t_POL);
     718       16359 :   y[1] = evalsigne(1) | evalvarn(v);
     719       16359 :   y++;
     720       16359 :   for (j=1; j<lx; j++) gel(y,j) = RgV_to_RgX(gel(x,j), w);
     721       16359 :   return normalizepol_lg(--y, lx+1);
     722             : }
     723             : 
     724             : /* matrix whose entries are given by the coeffs of the polynomial v in
     725             :  * two variables (considered as degree n-1 polynomials) */
     726             : GEN
     727         308 : RgXX_to_RgM(GEN v, long n)
     728             : {
     729         308 :   long j, N = lg(v)-1;
     730         308 :   GEN y = cgetg(N, t_MAT);
     731         308 :   for (j=1; j<N; j++) gel(y,j) = Rg_to_RgC(gel(v,j+1), n);
     732         308 :   return y;
     733             : }
     734             : 
     735             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     736             : GEN
     737       18918 : RgXY_swapspec(GEN x, long n, long w, long nx)
     738             : {
     739       18918 :   long j, ly = n+3;
     740       18918 :   GEN y = cgetg(ly, t_POL);
     741       18918 :   y[1] = evalsigne(1);
     742      243406 :   for (j=2; j<ly; j++)
     743             :   {
     744             :     long k;
     745      224488 :     GEN a = cgetg(nx+2,t_POL);
     746      224488 :     a[1] = evalsigne(1) | evalvarn(w);
     747     1152479 :     for (k=0; k<nx; k++)
     748             :     {
     749      927991 :       GEN xk = gel(x,k);
     750      927991 :       if (typ(xk)==t_POL)
     751      852734 :         gel(a,k+2) = j<lg(xk)? gel(xk,j): gen_0;
     752             :       else
     753       75257 :         gel(a,k+2) = j==2 ? xk: gen_0;
     754             :     }
     755      224488 :     gel(y,j) = normalizepol_lg(a, nx+2);
     756             :   }
     757       18918 :   return normalizepol_lg(y,ly);
     758             : }
     759             : 
     760             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     761             : GEN
     762         952 : RgXY_swap(GEN x, long n, long w)
     763             : {
     764         952 :   GEN z = RgXY_swapspec(x+2, n, w, lgpol(x));
     765         952 :   setvarn(z, varn(x)); return z;
     766             : }
     767             : 
     768             : long
     769         162 : RgXY_degreex(GEN b)
     770             : {
     771         162 :   long deg = -1, i;
     772         162 :   if (!signe(b)) return -1;
     773         976 :   for (i = 2; i < lg(b); ++i)
     774             :   {
     775         814 :     GEN bi = gel(b, i);
     776         814 :     if (typ(bi) == t_POL)
     777         792 :       deg = maxss(deg, degpol(bi));
     778             :   }
     779         162 :   return deg;
     780             : }
     781             : 
     782             : /* return (x % X^n). Shallow */
     783             : GEN
     784     2242403 : RgXn_red_shallow(GEN a, long n)
     785             : {
     786     2242403 :   long i, L = n+2, l = lg(a);
     787             :   GEN  b;
     788     2242403 :   if (L >= l) return a; /* deg(x) < n */
     789     1544790 :   b = cgetg(L, t_POL); b[1] = a[1];
     790     1544790 :   for (i=2; i<L; i++) gel(b,i) = gel(a,i);
     791     1544790 :   return normalizepol_lg(b,L);
     792             : }
     793             : 
     794             : GEN
     795         357 : RgXnV_red_shallow(GEN x, long n)
     796         357 : { pari_APPLY_type(t_VEC, RgXn_red_shallow(gel(x,i), n)) }
     797             : 
     798             : /* return (x * X^n). Shallow */
     799             : GEN
     800    60503672 : RgX_shift_shallow(GEN a, long n)
     801             : {
     802    60503672 :   long i, l = lg(a);
     803             :   GEN  b;
     804    60503672 :   if (l == 2 || !n) return a;
     805    43792761 :   l += n;
     806    43792761 :   if (n < 0)
     807             :   {
     808    38018844 :     if (l <= 2) return pol_0(varn(a));
     809    37995590 :     b = cgetg(l, t_POL); b[1] = a[1];
     810    37995591 :     a -= n;
     811    37995591 :     for (i=2; i<l; i++) gel(b,i) = gel(a,i);
     812             :   } else {
     813     5773917 :     b = cgetg(l, t_POL); b[1] = a[1];
     814     5773919 :     a -= n; n += 2;
     815     5773919 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     816     5773919 :     for (   ; i<l; i++) gel(b,i) = gel(a,i);
     817             :   }
     818    43769510 :   return b;
     819             : }
     820             : /* return (x * X^n). */
     821             : GEN
     822     3473034 : RgX_shift(GEN a, long n)
     823             : {
     824     3473034 :   long i, l = lg(a);
     825             :   GEN  b;
     826     3473034 :   if (l == 2 || !n) return RgX_copy(a);
     827     3472768 :   l += n;
     828     3472768 :   if (n < 0)
     829             :   {
     830         812 :     if (l <= 2) return pol_0(varn(a));
     831         770 :     b = cgetg(l, t_POL); b[1] = a[1];
     832         770 :     a -= n;
     833         770 :     for (i=2; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     834             :   } else {
     835     3471956 :     b = cgetg(l, t_POL); b[1] = a[1];
     836     3471956 :     a -= n; n += 2;
     837     3471956 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     838     3471956 :     for (   ; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     839             :   }
     840     3472726 :   return b;
     841             : }
     842             : 
     843             : GEN
     844      315273 : RgX_rotate_shallow(GEN P, long k, long p)
     845             : {
     846      315273 :   long i, l = lgpol(P);
     847             :   GEN r;
     848      315273 :   if (signe(P)==0)
     849         476 :     return pol_0(varn(P));
     850      314797 :   r = cgetg(p+2,t_POL); r[1] = P[1];
     851     2090746 :   for(i=0; i<p; i++)
     852             :   {
     853     1775949 :     long s = 2+(i+k)%p;
     854     1775949 :     gel(r,s) = i<l? gel(P,2+i): gen_0;
     855             :   }
     856      314797 :   return RgX_renormalize(r);
     857             : }
     858             : 
     859             : GEN
     860     2974324 : RgX_mulXn(GEN x, long d)
     861             : {
     862             :   pari_sp av;
     863             :   GEN z;
     864             :   long v;
     865     2974324 :   if (d >= 0) return RgX_shift(x, d);
     866     1463359 :   d = -d;
     867     1463359 :   v = RgX_val(x);
     868     1463359 :   if (v >= d) return RgX_shift(x, -d);
     869     1463345 :   av = avma;
     870     1463345 :   z = gred_rfrac_simple(RgX_shift_shallow(x, -v), pol_xn(d - v, varn(x)));
     871     1463345 :   return gerepileupto(av, z);
     872             : }
     873             : 
     874             : long
     875     2484334 : RgX_val(GEN x)
     876             : {
     877     2484334 :   long i, lx = lg(x);
     878     2484334 :   if (lx == 2) return LONG_MAX;
     879     2516555 :   for (i = 2; i < lx; i++)
     880     2516513 :     if (!isexactzero(gel(x,i))) break;
     881     2484292 :   if (i == lx) return LONG_MAX;/* possible with non-rational zeros */
     882     2484250 :   return i - 2;
     883             : }
     884             : long
     885    43956118 : RgX_valrem(GEN x, GEN *Z)
     886             : {
     887    43956118 :   long v, i, lx = lg(x);
     888    43956118 :   if (lx == 2) { *Z = pol_0(varn(x)); return LONG_MAX; }
     889    84089299 :   for (i = 2; i < lx; i++)
     890    84089299 :     if (!isexactzero(gel(x,i))) break;
     891             :   /* possible with non-rational zeros */
     892    43956118 :   if (i == lx) { *Z = pol_0(varn(x)); return LONG_MAX; }
     893    43956118 :   v = i - 2;
     894    43956118 :   *Z = RgX_shift_shallow(x, -v);
     895    43956118 :   return v;
     896             : }
     897             : long
     898       14528 : RgX_valrem_inexact(GEN x, GEN *Z)
     899             : {
     900             :   long v;
     901       14528 :   if (!signe(x)) { if (Z) *Z = pol_0(varn(x)); return LONG_MAX; }
     902       15018 :   for (v = 0;; v++)
     903       15515 :     if (!gequal0(gel(x,2+v))) break;
     904       14521 :   if (Z) *Z = RgX_shift_shallow(x, -v);
     905       14521 :   return v;
     906             : }
     907             : 
     908             : GEN
     909       89992 : RgXQC_red(GEN x, GEN T)
     910       89992 : { pari_APPLY_type(t_COL, grem(gel(x,i), T)) }
     911             : 
     912             : GEN
     913         945 : RgXQV_red(GEN x, GEN T)
     914         945 : { pari_APPLY_type(t_VEC, grem(gel(x,i), T)) }
     915             : 
     916             : GEN
     917       13090 : RgXQM_red(GEN x, GEN T)
     918       13090 : { pari_APPLY_same(RgXQC_red(gel(x,i), T)) }
     919             : 
     920             : GEN
     921           0 : RgXQM_mul(GEN P, GEN Q, GEN T)
     922             : {
     923           0 :   return RgXQM_red(RgM_mul(P, Q), T);
     924             : }
     925             : 
     926             : GEN
     927      379676 : RgXQX_red(GEN P, GEN T)
     928             : {
     929      379676 :   long i, l = lg(P);
     930      379676 :   GEN Q = cgetg(l, t_POL);
     931      379676 :   Q[1] = P[1];
     932      379676 :   for (i=2; i<l; i++) gel(Q,i) = grem(gel(P,i), T);
     933      379676 :   return normalizepol_lg(Q, l);
     934             : }
     935             : 
     936             : GEN
     937      239334 : RgX_deriv(GEN x)
     938             : {
     939      239334 :   long i,lx = lg(x)-1;
     940             :   GEN y;
     941             : 
     942      239334 :   if (lx<3) return pol_0(varn(x));
     943      211278 :   y = cgetg(lx,t_POL); gel(y,2) = gcopy(gel(x,3));
     944      211278 :   for (i=3; i<lx ; i++) gel(y,i) = gmulsg(i-1,gel(x,i+1));
     945      211278 :   y[1] = x[1]; return normalizepol_lg(y,i);
     946             : }
     947             : 
     948             : GEN
     949      477156 : RgX_recipspec_shallow(GEN x, long l, long n)
     950             : {
     951             :   long i;
     952      477156 :   GEN z = cgetg(n+2,t_POL);
     953      477162 :   z[1] = 0; z += 2;
     954    14476367 :   for(i=0; i<l; i++)
     955    13999205 :     gel(z,n-i-1) = gel(x,i);
     956      604201 :   for(   ; i<n; i++)
     957      127039 :     gel(z, n-i-1) = gen_0;
     958      477162 :   return normalizepol_lg(z-2,n+2);
     959             : }
     960             : 
     961             : GEN
     962       40979 : RgXn_recip_shallow(GEN P, long n)
     963             : {
     964       40979 :   GEN Q = RgX_recipspec_shallow(P+2, lgpol(P), n);
     965       40979 :   setvarn(Q, varn(P));
     966       40979 :   return Q;
     967             : }
     968             : 
     969             : /* return coefficients s.t x = x_0 X^n + ... + x_n */
     970             : GEN
     971        2100 : RgX_recip(GEN x)
     972             : {
     973             :   long lx, i, j;
     974        2100 :   GEN y = cgetg_copy(x, &lx);
     975        2100 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gcopy(gel(x,j));
     976        2100 :   return normalizepol_lg(y,lx);
     977             : }
     978             : /* shallow version */
     979             : GEN
     980      376846 : RgX_recip_shallow(GEN x)
     981             : {
     982             :   long lx, i, j;
     983      376846 :   GEN y = cgetg_copy(x, &lx);
     984      376885 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
     985      376885 :   return y;
     986             : }
     987             : /*******************************************************************/
     988             : /*                                                                 */
     989             : /*                      ADDITION / SUBTRACTION                     */
     990             : /*                                                                 */
     991             : /*******************************************************************/
     992             : /* same variable */
     993             : GEN
     994    18226294 : RgX_add(GEN x, GEN y)
     995             : {
     996    18226294 :   long i, lx = lg(x), ly = lg(y);
     997             :   GEN z;
     998    18226294 :   if (ly <= lx) {
     999    16052349 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1000    16052349 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1001    16052349 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1002    16052349 :     z = normalizepol_lg(z, lx);
    1003             :   } else {
    1004     2173945 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1005     2173945 :     for (i=2; i < lx; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1006     2173945 :     for (   ; i < ly; i++) gel(z,i) = gcopy(gel(y,i));
    1007     2173945 :     z = normalizepol_lg(z, ly);
    1008             :   }
    1009    18226294 :   return z;
    1010             : }
    1011             : GEN
    1012    12825979 : RgX_sub(GEN x, GEN y)
    1013             : {
    1014    12825979 :   long i, lx = lg(x), ly = lg(y);
    1015             :   GEN z;
    1016    12825979 :   if (ly <= lx) {
    1017     9863113 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1018     9863113 :     for (i=2; i < ly; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1019     9863113 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1020     9863113 :     z = normalizepol_lg(z, lx);
    1021             :   } else {
    1022     2962866 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1023     2962866 :     for (i=2; i < lx; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1024     2962866 :     for (   ; i < ly; i++) gel(z,i) = gneg(gel(y,i));
    1025     2962866 :     z = normalizepol_lg(z, ly);
    1026             :   }
    1027    12825979 :   return z;
    1028             : }
    1029             : GEN
    1030     1098393 : RgX_neg(GEN x)
    1031             : {
    1032     1098393 :   long i, lx = lg(x);
    1033     1098393 :   GEN y = cgetg(lx, t_POL); y[1] = x[1];
    1034     1098393 :   for (i=2; i<lx; i++) gel(y,i) = gneg(gel(x,i));
    1035     1098393 :   return y;
    1036             : }
    1037             : 
    1038             : GEN
    1039    12614870 : RgX_Rg_add(GEN y, GEN x)
    1040             : {
    1041             :   GEN z;
    1042    12614870 :   long lz = lg(y), i;
    1043    12614870 :   if (lz == 2) return scalarpol(x,varn(y));
    1044     9659614 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1045     9659614 :   gel(z,2) = gadd(gel(y,2),x);
    1046     9659614 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1047             :   /* probably useless unless lz = 3, but cannot be skipped if y is
    1048             :    * an inexact 0 */
    1049     9659614 :   return normalizepol_lg(z,lz);
    1050             : }
    1051             : GEN
    1052        2548 : RgX_Rg_add_shallow(GEN y, GEN x)
    1053             : {
    1054             :   GEN z;
    1055        2548 :   long lz = lg(y), i;
    1056        2548 :   if (lz == 2) return scalarpol(x,varn(y));
    1057        2548 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1058        2548 :   gel(z,2) = gadd(gel(y,2),x);
    1059        2548 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
    1060        2548 :   return z = normalizepol_lg(z,lz);
    1061             : }
    1062             : GEN
    1063       31095 : RgX_Rg_sub(GEN y, GEN x)
    1064             : {
    1065             :   GEN z;
    1066       31095 :   long lz = lg(y), i;
    1067       31095 :   if (lz == 2)
    1068             :   { /* scalarpol(gneg(x),varn(y)) optimized */
    1069        2576 :     long v = varn(y);
    1070        2576 :     if (isrationalzero(x)) return pol_0(v);
    1071          14 :     z = cgetg(3,t_POL);
    1072          28 :     z[1] = gequal0(x)? evalvarn(v)
    1073          14 :                    : evalvarn(v) | evalsigne(1);
    1074          14 :     gel(z,2) = gneg(x); return z;
    1075             :   }
    1076       28519 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1077       28519 :   gel(z,2) = gsub(gel(y,2),x);
    1078       28519 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1079       28519 :   return z = normalizepol_lg(z,lz);
    1080             : }
    1081             : GEN
    1082      469448 : Rg_RgX_sub(GEN x, GEN y)
    1083             : {
    1084             :   GEN z;
    1085      469448 :   long lz = lg(y), i;
    1086      469448 :   if (lz == 2) return scalarpol(x,varn(y));
    1087      468342 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1088      468342 :   gel(z,2) = gsub(x, gel(y,2));
    1089      468342 :   for(i=3; i<lz; i++) gel(z,i) = gneg(gel(y,i));
    1090      468342 :   return z = normalizepol_lg(z,lz);
    1091             : }
    1092             : /*******************************************************************/
    1093             : /*                                                                 */
    1094             : /*                  KARATSUBA MULTIPLICATION                       */
    1095             : /*                                                                 */
    1096             : /*******************************************************************/
    1097             : #if 0
    1098             : /* to debug Karatsuba-like routines */
    1099             : GEN
    1100             : zx_debug_spec(GEN x, long nx)
    1101             : {
    1102             :   GEN z = cgetg(nx+2,t_POL);
    1103             :   long i;
    1104             :   for (i=0; i<nx; i++) gel(z,i+2) = stoi(x[i]);
    1105             :   z[1] = evalsigne(1); return z;
    1106             : }
    1107             : 
    1108             : GEN
    1109             : RgX_debug_spec(GEN x, long nx)
    1110             : {
    1111             :   GEN z = cgetg(nx+2,t_POL);
    1112             :   long i;
    1113             :   for (i=0; i<nx; i++) z[i+2] = x[i];
    1114             :   z[1] = evalsigne(1); return z;
    1115             : }
    1116             : #endif
    1117             : 
    1118             : /* generic multiplication */
    1119             : GEN
    1120     3194712 : RgX_addspec_shallow(GEN x, GEN y, long nx, long ny)
    1121             : {
    1122             :   GEN z, t;
    1123             :   long i;
    1124     3194712 :   if (nx == ny) {
    1125      744124 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1126      744124 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1127      744124 :     return normalizepol_lg(z, nx+2);
    1128             :   }
    1129     2450588 :   if (ny < nx) {
    1130     2298640 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1131     2298640 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1132     2298640 :     for (   ; i < nx; i++) gel(t,i) = gel(x,i);
    1133     2298640 :     return normalizepol_lg(z, nx+2);
    1134             :   } else {
    1135      151948 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1136      151948 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1137      151948 :     for (   ; i < ny; i++) gel(t,i) = gel(y,i);
    1138      151948 :     return normalizepol_lg(z, ny+2);
    1139             :   }
    1140             : }
    1141             : GEN
    1142      203113 : RgX_addspec(GEN x, GEN y, long nx, long ny)
    1143             : {
    1144             :   GEN z, t;
    1145             :   long i;
    1146      203113 :   if (nx == ny) {
    1147       10423 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1148       10423 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1149       10423 :     return normalizepol_lg(z, nx+2);
    1150             :   }
    1151      192690 :   if (ny < nx) {
    1152      190989 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1153      190989 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1154      190989 :     for (   ; i < nx; i++) gel(t,i) = gcopy(gel(x,i));
    1155      190989 :     return normalizepol_lg(z, nx+2);
    1156             :   } else {
    1157        1701 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1158        1701 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1159        1701 :     for (   ; i < ny; i++) gel(t,i) = gcopy(gel(y,i));
    1160        1701 :     return normalizepol_lg(z, ny+2);
    1161             :   }
    1162             : }
    1163             : 
    1164             : /* Return the vector of coefficients of x, where we replace rational 0s by NULL
    1165             :  * [ to speed up basic operation s += x[i]*y[j] ]. We create a proper
    1166             :  * t_VECSMALL, to hold this, which can be left on stack: gerepile
    1167             :  * will not crash on it. The returned vector itself is not a proper GEN,
    1168             :  * we access the coefficients as x[i], i = 0..deg(x) */
    1169             : static GEN
    1170    20960034 : RgXspec_kill0(GEN x, long lx)
    1171             : {
    1172    20960034 :   GEN z = cgetg(lx+1, t_VECSMALL) + 1; /* inhibit gerepile-wise */
    1173             :   long i;
    1174    84200794 :   for (i=0; i <lx; i++)
    1175             :   {
    1176    63240760 :     GEN c = gel(x,i);
    1177    63240760 :     z[i] = (long)(isrationalzero(c)? NULL: c);
    1178             :   }
    1179    20960034 :   return z;
    1180             : }
    1181             : 
    1182             : INLINE GEN
    1183    46439569 : RgX_mulspec_basecase_limb(GEN x, GEN y, long a, long b)
    1184             : {
    1185    46439569 :   pari_sp av = avma;
    1186    46439569 :   GEN s = NULL;
    1187             :   long i;
    1188             : 
    1189   208864899 :   for (i=a; i<b; i++)
    1190   162425330 :     if (gel(y,i) && gel(x,-i))
    1191             :     {
    1192    86025406 :       GEN t = gmul(gel(y,i), gel(x,-i));
    1193    86025406 :       s = s? gadd(s, t): t;
    1194             :     }
    1195    46439569 :   return s? gerepileupto(av, s): gen_0;
    1196             : }
    1197             : 
    1198             : /* assume nx >= ny > 0, return x * y * t^v */
    1199             : static GEN
    1200     8255721 : RgX_mulspec_basecase(GEN x, GEN y, long nx, long ny, long v)
    1201             : {
    1202             :   long i, lz, nz;
    1203             :   GEN z;
    1204             : 
    1205     8255721 :   x = RgXspec_kill0(x,nx);
    1206     8255721 :   y = RgXspec_kill0(y,ny);
    1207     8255721 :   lz = nx + ny + 1; nz = lz-2;
    1208     8255721 :   lz += v;
    1209     8255721 :   z = cgetg(lz, t_POL) + 2; /* x:y:z [i] = term of degree i */
    1210     8255721 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1211     8255721 :   for (i=0; i<ny; i++)gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0, i+1);
    1212     8255721 :   for (  ; i<nx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ny);
    1213     8255721 :   for (  ; i<nz; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-nx+1,ny);
    1214     8255721 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1215             : }
    1216             : 
    1217             : /* return (x * X^d) + y. Assume d > 0 */
    1218             : GEN
    1219     2322412 : RgX_addmulXn_shallow(GEN x0, GEN y0, long d)
    1220             : {
    1221             :   GEN x, y, xd, yd, zd;
    1222             :   long a, lz, nx, ny;
    1223             : 
    1224     2322412 :   if (!signe(x0)) return y0;
    1225     2298463 :   ny = lgpol(y0);
    1226     2298463 :   nx = lgpol(x0);
    1227     2298463 :   zd = (GEN)avma;
    1228     2298463 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1229     2298463 :   if (a <= 0)
    1230             :   {
    1231      184509 :     lz = nx+d+2;
    1232      184509 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1233      184509 :     while (xd > x) gel(--zd,0) = gel(--xd,0);
    1234      184509 :     x = zd + a;
    1235      184509 :     while (zd > x) gel(--zd,0) = gen_0;
    1236             :   }
    1237             :   else
    1238             :   {
    1239     2113954 :     xd = new_chunk(d); yd = y+d;
    1240     2113954 :     x = RgX_addspec_shallow(x,yd, nx,a);
    1241     2113954 :     lz = (a>nx)? ny+2: lg(x)+d;
    1242     2113954 :     x += 2; while (xd > x) *--zd = *--xd;
    1243             :   }
    1244     2298463 :   while (yd > y) *--zd = *--yd;
    1245     2298463 :   *--zd = x0[1];
    1246     2298463 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1247             : }
    1248             : GEN
    1249      495073 : RgX_addmulXn(GEN x0, GEN y0, long d)
    1250             : {
    1251             :   GEN x, y, xd, yd, zd;
    1252             :   long a, lz, nx, ny;
    1253             : 
    1254      495073 :   if (!signe(x0)) return RgX_copy(y0);
    1255      494496 :   nx = lgpol(x0);
    1256      494496 :   ny = lgpol(y0);
    1257      494496 :   zd = (GEN)avma;
    1258      494496 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1259      494496 :   if (a <= 0)
    1260             :   {
    1261      291383 :     lz = nx+d+2;
    1262      291383 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1263      291383 :     while (xd > x) gel(--zd,0) = gcopy(gel(--xd,0));
    1264      291383 :     x = zd + a;
    1265      291383 :     while (zd > x) gel(--zd,0) = gen_0;
    1266             :   }
    1267             :   else
    1268             :   {
    1269      203113 :     xd = new_chunk(d); yd = y+d;
    1270      203113 :     x = RgX_addspec(x,yd, nx,a);
    1271      203113 :     lz = (a>nx)? ny+2: lg(x)+d;
    1272      203113 :     x += 2; while (xd > x) *--zd = *--xd;
    1273             :   }
    1274      494496 :   while (yd > y) gel(--zd,0) = gcopy(gel(--yd,0));
    1275      494496 :   *--zd = x0[1];
    1276      494496 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1277             : }
    1278             : 
    1279             : /* return x * y mod t^n */
    1280             : static GEN
    1281     2221375 : RgXn_mul_basecase(GEN x, GEN y, long n)
    1282             : {
    1283     2221375 :   long i, lz = n+2, lx = lgpol(x), ly = lgpol(y);
    1284             :   GEN z;
    1285     2221375 :   if (lx < 0) return pol_0(varn(x));
    1286     2221375 :   if (ly < 0) return pol_0(varn(x));
    1287     2221375 :   z = cgetg(lz, t_POL) + 2;
    1288     2221375 :   x+=2; if (lx > n) lx = n;
    1289     2221375 :   y+=2; if (ly > n) ly = n;
    1290     2221375 :   z[-1] = x[-1];
    1291     2221375 :   if (ly > lx) { swap(x,y); lswap(lx,ly); }
    1292     2221375 :   x = RgXspec_kill0(x, lx);
    1293     2221375 :   y = RgXspec_kill0(y, ly);
    1294             :   /* x:y:z [i] = term of degree i */
    1295     2221375 :   for (i=0;i<ly; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,i+1);
    1296     2221375 :   for (  ; i<lx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ly);
    1297     2221375 :   for (  ; i<n; i++)  gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-lx+1,ly);
    1298     2221375 :   return normalizepol_lg(z - 2, lz);
    1299             : }
    1300             : /* Mulders / Karatsuba product f*g mod t^n (Hanrot-Zimmermann variant) */
    1301             : static GEN
    1302     2306796 : RgXn_mul2(GEN f, GEN g, long n)
    1303             : {
    1304     2306796 :   pari_sp av = avma;
    1305             :   GEN fe,fo, ge,go, l,h,m;
    1306             :   long n0, n1;
    1307     2306796 :   if (degpol(f) + degpol(g) < n) return RgX_mul(f,g);
    1308     2247289 :   if (n < 80) return RgXn_mul_basecase(f,g,n);
    1309       25914 :   n0 = n>>1; n1 = n-n0;
    1310       25914 :   RgX_even_odd(f, &fe, &fo);
    1311       25914 :   RgX_even_odd(g, &ge, &go);
    1312       25914 :   l = RgXn_mul(fe,ge,n1);
    1313       25914 :   h = RgXn_mul(fo,go,n0);
    1314       25914 :   m = RgX_sub(RgXn_mul(RgX_add(fe,fo),RgX_add(ge,go),n0), RgX_add(l,h));
    1315             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1316             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1317       25914 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1318             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1319       25914 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1320       25914 :   m = RgX_inflate(m,2);
    1321             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1322       25914 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1323       25914 :   h = RgX_inflate(h,2);
    1324       25914 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1325       25914 :   return gerepileupto(av, h);
    1326             : }
    1327             : /* (f*g) \/ x^n */
    1328             : static GEN
    1329        2470 : RgX_mulhigh_i2(GEN f, GEN g, long n)
    1330             : {
    1331        2470 :   long d = degpol(f)+degpol(g) + 1 - n;
    1332             :   GEN h;
    1333        2470 :   if (d <= 2) return RgX_shift_shallow(RgX_mul(f,g), -n);
    1334        1669 :   h = RgX_recip_shallow(RgXn_mul(RgX_recip_shallow(f),
    1335             :                                  RgX_recip_shallow(g), d));
    1336        1669 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1337             : }
    1338             : 
    1339             : /* (f*g) \/ x^n */
    1340             : static GEN
    1341           0 : RgX_sqrhigh_i2(GEN f, long n)
    1342             : {
    1343           0 :   long d = 2*degpol(f)+ 1 - n;
    1344             :   GEN h;
    1345           0 :   if (d <= 2) return RgX_shift_shallow(RgX_sqr(f), -n);
    1346           0 :   h = RgX_recip_shallow(RgXn_sqr(RgX_recip_shallow(f), d));
    1347           0 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1348             : }
    1349             : 
    1350             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
    1351             :  * b+2 were sent instead. na, nb = number of terms of a, b.
    1352             :  * Only c, c0, c1, c2 are genuine GEN.
    1353             :  */
    1354             : GEN
    1355     8824994 : RgX_mulspec(GEN a, GEN b, long na, long nb)
    1356             : {
    1357             :   GEN a0, c, c0;
    1358     8824994 :   long n0, n0a, i, v = 0;
    1359             :   pari_sp av;
    1360             : 
    1361     8824994 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v++; }
    1362     8824994 :   while (nb && isrationalzero(gel(b,0))) { b++; nb--; v++; }
    1363     8824994 :   if (na < nb) swapspec(a,b, na,nb);
    1364     8824994 :   if (!nb) return pol_0(0);
    1365             : 
    1366     8724841 :   if (nb < RgX_MUL_LIMIT) return RgX_mulspec_basecase(a,b,na,nb, v);
    1367      469120 :   RgX_shift_inplace_init(v);
    1368      469120 :   i = (na>>1); n0 = na-i; na = i;
    1369      469120 :   av = avma; a0 = a+n0; n0a = n0;
    1370      469120 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1371             : 
    1372      469120 :   if (nb > n0)
    1373             :   {
    1374             :     GEN b0,c1,c2;
    1375             :     long n0b;
    1376             : 
    1377      467796 :     nb -= n0; b0 = b+n0; n0b = n0;
    1378      467796 :     while (n0b && isrationalzero(gel(b,n0b-1))) n0b--;
    1379      467796 :     c = RgX_mulspec(a,b,n0a,n0b);
    1380      467796 :     c0 = RgX_mulspec(a0,b0, na,nb);
    1381             : 
    1382      467796 :     c2 = RgX_addspec_shallow(a0,a, na,n0a);
    1383      467796 :     c1 = RgX_addspec_shallow(b0,b, nb,n0b);
    1384             : 
    1385      467796 :     c1 = RgX_mulspec(c1+2,c2+2, lgpol(c1),lgpol(c2));
    1386      467796 :     c2 = RgX_sub(c1, RgX_add(c0,c));
    1387      467796 :     c0 = RgX_addmulXn_shallow(c0, c2, n0);
    1388             :   }
    1389             :   else
    1390             :   {
    1391        1324 :     c = RgX_mulspec(a,b,n0a,nb);
    1392        1324 :     c0 = RgX_mulspec(a0,b,na,nb);
    1393             :   }
    1394      469120 :   c0 = RgX_addmulXn(c0,c,n0);
    1395      469120 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1396             : }
    1397             : 
    1398             : INLINE GEN
    1399       33134 : RgX_sqrspec_basecase_limb(GEN x, long a, long i)
    1400             : {
    1401       33134 :   pari_sp av = avma;
    1402       33134 :   GEN s = NULL;
    1403       33134 :   long j, l = (i+1)>>1;
    1404       98812 :   for (j=a; j<l; j++)
    1405             :   {
    1406       65678 :     GEN xj = gel(x,j), xx = gel(x,i-j);
    1407       65678 :     if (xj && xx)
    1408             :     {
    1409       58356 :       GEN t = gmul(xj, xx);
    1410       58356 :       s = s? gadd(s, t): t;
    1411             :     }
    1412             :   }
    1413       33134 :   if (s) s = gshift(s,1);
    1414       33134 :   if ((i&1) == 0)
    1415             :   {
    1416       19488 :     GEN t = gel(x, i>>1);
    1417       19488 :     if (t) {
    1418       17171 :       t = gsqr(t);
    1419       17171 :       s = s? gadd(s, t): t;
    1420             :     }
    1421             :   }
    1422       33134 :   return s? gerepileupto(av,s): gen_0;
    1423             : }
    1424             : static GEN
    1425        5842 : RgX_sqrspec_basecase(GEN x, long nx, long v)
    1426             : {
    1427             :   long i, lz, nz;
    1428             :   GEN z;
    1429             : 
    1430        5842 :   if (!nx) return pol_0(0);
    1431        5842 :   x = RgXspec_kill0(x,nx);
    1432        5842 :   lz = (nx << 1) + 1, nz = lz-2;
    1433        5842 :   lz += v;
    1434        5842 :   z = cgetg(lz,t_POL) + 2;
    1435        5842 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1436        5842 :   for (i=0; i<nx; i++)gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1437        5842 :   for (  ; i<nz; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, i-nx+1, i);
    1438        5842 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1439             : }
    1440             : /* return x^2 mod t^n */
    1441             : static GEN
    1442           0 : RgXn_sqr_basecase(GEN x, long n)
    1443             : {
    1444           0 :   long i, lz = n+2, lx = lgpol(x);
    1445             :   GEN z;
    1446           0 :   if (lx < 0) return pol_0(varn(x));
    1447           0 :   z = cgetg(lz, t_POL);
    1448           0 :   z[1] = x[1];
    1449           0 :   x+=2; if (lx > n) lx = n;
    1450           0 :   x = RgXspec_kill0(x,lx);
    1451           0 :   z+=2;/* x:z [i] = term of degree i */
    1452           0 :   for (i=0;i<lx; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1453           0 :   for (  ; i<n; i++)  gel(z,i) = RgX_sqrspec_basecase_limb(x, i-lx+1, i);
    1454           0 :   z -= 2; return normalizepol_lg(z, lz);
    1455             : }
    1456             : /* Mulders / Karatsuba product f^2 mod t^n (Hanrot-Zimmermann variant) */
    1457             : static GEN
    1458         238 : RgXn_sqr2(GEN f, long n)
    1459             : {
    1460         238 :   pari_sp av = avma;
    1461             :   GEN fe,fo, l,h,m;
    1462             :   long n0, n1;
    1463         238 :   if (2*degpol(f) < n) return RgX_sqr_i(f);
    1464           0 :   if (n < 80) return RgXn_sqr_basecase(f,n);
    1465           0 :   n0 = n>>1; n1 = n-n0;
    1466           0 :   RgX_even_odd(f, &fe, &fo);
    1467           0 :   l = RgXn_sqr(fe,n1);
    1468           0 :   h = RgXn_sqr(fo,n0);
    1469           0 :   m = RgX_sub(RgXn_sqr(RgX_add(fe,fo),n0), RgX_add(l,h));
    1470             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1471             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1472           0 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1473             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1474           0 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1475           0 :   m = RgX_inflate(m,2);
    1476             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1477           0 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1478           0 :   h = RgX_inflate(h,2);
    1479           0 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1480           0 :   return gerepileupto(av, h);
    1481             : }
    1482             : GEN
    1483        5881 : RgX_sqrspec(GEN a, long na)
    1484             : {
    1485             :   GEN a0, c, c0, c1;
    1486        5881 :   long n0, n0a, i, v = 0;
    1487             :   pari_sp av;
    1488             : 
    1489        5881 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v += 2; }
    1490        5881 :   if (na<RgX_SQR_LIMIT) return RgX_sqrspec_basecase(a, na, v);
    1491          39 :   RgX_shift_inplace_init(v);
    1492          39 :   i = (na>>1); n0 = na-i; na = i;
    1493          39 :   av = avma; a0 = a+n0; n0a = n0;
    1494          39 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1495             : 
    1496          39 :   c = RgX_sqrspec(a,n0a);
    1497          39 :   c0 = RgX_sqrspec(a0,na);
    1498          39 :   c1 = gmul2n(RgX_mulspec(a0,a, na,n0a), 1);
    1499          39 :   c0 = RgX_addmulXn_shallow(c0,c1, n0);
    1500          39 :   c0 = RgX_addmulXn(c0,c,n0);
    1501          39 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1502             : }
    1503             : 
    1504             : /* (X^a + A)(X^b + B) - X^(a+b), where deg A < a, deg B < b */
    1505             : GEN
    1506      689928 : RgX_mul_normalized(GEN A, long a, GEN B, long b)
    1507             : {
    1508      689928 :   GEN z = RgX_mul(A, B);
    1509      689928 :   if (a < b)
    1510        5446 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(A, B, b-a), z, a);
    1511      684482 :   else if (a > b)
    1512      432426 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(B, A, a-b), z, b);
    1513             :   else
    1514      252056 :     z = RgX_addmulXn_shallow(RgX_add(A, B), z, a);
    1515      689928 :   return z;
    1516             : }
    1517             : 
    1518             : GEN
    1519     7418919 : RgX_mul_i(GEN x, GEN y)
    1520             : {
    1521     7418919 :   GEN z = RgX_mulspec(x+2, y+2, lgpol(x), lgpol(y));
    1522     7418919 :   setvarn(z, varn(x)); return z;
    1523             : }
    1524             : 
    1525             : GEN
    1526        5803 : RgX_sqr_i(GEN x)
    1527             : {
    1528        5803 :   GEN z = RgX_sqrspec(x+2, lgpol(x));
    1529        5803 :   setvarn(z,varn(x)); return z;
    1530             : }
    1531             : 
    1532             : /*******************************************************************/
    1533             : /*                                                                 */
    1534             : /*                               DIVISION                          */
    1535             : /*                                                                 */
    1536             : /*******************************************************************/
    1537             : GEN
    1538      438856 : RgX_Rg_divexact(GEN x, GEN y) {
    1539             :   long i, lx;
    1540             :   GEN z;
    1541      438856 :   if (typ(y) == t_INT && is_pm1(y))
    1542       28124 :     return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1543      410732 :   z = cgetg_copy(x, &lx); z[1] = x[1];
    1544      410732 :   for (i=2; i<lx; i++) gel(z,i) = gdivexact(gel(x,i),y);
    1545      410732 :   return z;
    1546             : }
    1547             : GEN
    1548    23885883 : RgX_Rg_div(GEN x, GEN y) {
    1549             :   long i, lx;
    1550    23885883 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1551    23885883 :   for (i=2; i<lx; i++) gel(z,i) = gdiv(gel(x,i),y);
    1552    23885883 :   return normalizepol_lg(z, lx);
    1553             : }
    1554             : GEN
    1555       12698 : RgX_normalize(GEN x)
    1556             : {
    1557       12698 :   GEN z, d = NULL;
    1558       12698 :   long i, n = lg(x)-1;
    1559       12698 :   for (i = n; i > 1; i--) { d = gel(x,i); if (!gequal0(d)) break; }
    1560       12698 :   if (i == 1) return pol_0(varn(x));
    1561       12698 :   if (i == n && isint1(d)) return x;
    1562        5621 :   n = i; z = cgetg(n+1, t_POL); z[1] = x[1];
    1563        5621 :   for (i=2; i<n; i++) gel(z,i) = gdiv(gel(x,i),d);
    1564        5621 :   gel(z,n) = Rg_get_1(d); return z;
    1565             : }
    1566             : GEN
    1567        2597 : RgX_divs(GEN x, long y) {
    1568             :   long i, lx;
    1569        2597 :   GEN z = cgetg_copy(x, &lx); z[1] = x[1];
    1570        2597 :   for (i=2; i<lx; i++) gel(z,i) = gdivgs(gel(x,i),y);
    1571        2597 :   return normalizepol_lg(z, lx);
    1572             : }
    1573             : GEN
    1574       41228 : RgX_div_by_X_x(GEN a, GEN x, GEN *r)
    1575             : {
    1576       41228 :   long l = lg(a), i;
    1577       41228 :   GEN a0, z0, z = cgetg(l-1, t_POL);
    1578       41228 :   z[1] = a[1];
    1579       41228 :   a0 = a + l-1;
    1580       41228 :   z0 = z + l-2; *z0 = *a0--;
    1581      790828 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1582             :   {
    1583      749600 :     GEN t = gadd(gel(a0--,0), gmul(x, gel(z0--,0)));
    1584      749600 :     gel(z0,0) = t;
    1585             :   }
    1586       41228 :   if (r) *r = gadd(gel(a0,0), gmul(x, gel(z0,0)));
    1587       41228 :   return z;
    1588             : }
    1589             : /* Polynomial division x / y:
    1590             :  *   if pr = ONLY_REM return remainder, otherwise return quotient
    1591             :  *   if pr = ONLY_DIVIDES return quotient if division is exact, else NULL
    1592             :  *   if pr != NULL set *pr to remainder, as the last object on stack */
    1593             : /* assume, typ(x) = typ(y) = t_POL, same variable */
    1594             : static GEN
    1595    12214775 : RgX_divrem_i(GEN x, GEN y, GEN *pr)
    1596             : {
    1597             :   pari_sp avy, av, av1;
    1598             :   long dx,dy,dz,i,j,sx,lr;
    1599             :   GEN z,p1,p2,rem,y_lead,mod,p;
    1600             :   GEN (*f)(GEN,GEN);
    1601             : 
    1602    12214775 :   if (!signe(y)) pari_err_INV("RgX_divrem",y);
    1603             : 
    1604    12214775 :   dy = degpol(y);
    1605    12214775 :   y_lead = gel(y,dy+2);
    1606    12214775 :   if (gequal0(y_lead)) /* normalize denominator if leading term is 0 */
    1607             :   {
    1608           0 :     pari_warn(warner,"normalizing a polynomial with 0 leading term");
    1609           0 :     for (dy--; dy>=0; dy--)
    1610             :     {
    1611           0 :       y_lead = gel(y,dy+2);
    1612           0 :       if (!gequal0(y_lead)) break;
    1613             :     }
    1614             :   }
    1615    12214775 :   if (!dy) /* y is constant */
    1616             :   {
    1617       52665 :     if (pr == ONLY_REM) return pol_0(varn(x));
    1618       52665 :     z = RgX_Rg_div(x, y_lead);
    1619       52665 :     if (pr == ONLY_DIVIDES) return z;
    1620       50915 :     if (pr) *pr = pol_0(varn(x));
    1621       50915 :     return z;
    1622             :   }
    1623    12162110 :   dx = degpol(x);
    1624    12162110 :   if (dx < dy)
    1625             :   {
    1626      637784 :     if (pr == ONLY_REM) return RgX_copy(x);
    1627      371691 :     if (pr == ONLY_DIVIDES) return signe(x)? NULL: pol_0(varn(x));
    1628      371670 :     z = pol_0(varn(x));
    1629      371670 :     if (pr) *pr = RgX_copy(x);
    1630      371670 :     return z;
    1631             :   }
    1632             : 
    1633             :   /* x,y in R[X], y non constant */
    1634    11524326 :   av = avma;
    1635    11524326 :   p = NULL;
    1636    11524326 :   if (RgX_is_FpX(x, &p) && RgX_is_FpX(y, &p) && p)
    1637             :   {
    1638      225904 :     z = FpX_divrem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p, pr);
    1639      225904 :     if (!z) return gc_NULL(av);
    1640      225904 :     z = FpX_to_mod(z, p);
    1641      225904 :     if (!pr || pr == ONLY_REM || pr == ONLY_DIVIDES)
    1642      121681 :       return gerepileupto(av, z);
    1643      104223 :     *pr = FpX_to_mod(*pr, p);
    1644      104223 :     gerepileall(av, 2, pr, &z);
    1645      104223 :     return z;
    1646             :   }
    1647    11298422 :   switch(typ(y_lead))
    1648             :   {
    1649             :     case t_REAL:
    1650           0 :       y_lead = ginv(y_lead);
    1651           0 :       f = gmul; mod = NULL;
    1652           0 :       break;
    1653             :     case t_INTMOD:
    1654        2001 :     case t_POLMOD: y_lead = ginv(y_lead);
    1655        2001 :       f = gmul; mod = gmodulo(gen_1, gel(y_lead,1));
    1656        2001 :       break;
    1657    11296421 :     default: if (gequal1(y_lead)) y_lead = NULL;
    1658    11296421 :       f = gdiv; mod = NULL;
    1659             :   }
    1660             : 
    1661    11298422 :   if (y_lead == NULL)
    1662     9520009 :     p2 = gel(x,dx+2);
    1663             :   else {
    1664             :     for(;;) {
    1665     1778413 :       p2 = f(gel(x,dx+2),y_lead);
    1666     1778413 :       p2 = simplify_shallow(p2);
    1667     1778413 :       if (!isexactzero(p2) || (--dx < 0)) break;
    1668             :     }
    1669     1778413 :     if (dx < dy) /* leading coeff of x was in fact zero */
    1670             :     {
    1671           0 :       if (pr == ONLY_DIVIDES) {
    1672           0 :         set_avma(av);
    1673           0 :         return (dx < 0)? pol_0(varn(x)) : NULL;
    1674             :       }
    1675           0 :       if (pr == ONLY_REM)
    1676             :       {
    1677           0 :         if (dx < 0)
    1678           0 :           return gerepilecopy(av, scalarpol(p2, varn(x)));
    1679             :         else
    1680             :         {
    1681             :           GEN t;
    1682           0 :           set_avma(av);
    1683           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1684           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1685           0 :           return t;
    1686             :         }
    1687             :       }
    1688           0 :       if (pr) /* cf ONLY_REM above */
    1689             :       {
    1690           0 :         if (dx < 0)
    1691             :         {
    1692           0 :           p2 = gclone(p2);
    1693           0 :           set_avma(av);
    1694           0 :           z = pol_0(varn(x));
    1695           0 :           x = scalarpol(p2, varn(x));
    1696           0 :           gunclone(p2);
    1697             :         }
    1698             :         else
    1699             :         {
    1700             :           GEN t;
    1701           0 :           set_avma(av);
    1702           0 :           z = pol_0(varn(x));
    1703           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1704           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1705           0 :           x = t;
    1706             :         }
    1707           0 :         *pr = x;
    1708             :       }
    1709             :       else
    1710             :       {
    1711           0 :         set_avma(av);
    1712           0 :         z = pol_0(varn(x));
    1713             :       }
    1714           0 :       return z;
    1715             :     }
    1716             :   }
    1717             :   /* dx >= dy */
    1718    11298422 :   avy = avma;
    1719    11298422 :   dz = dx-dy;
    1720    11298422 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1721    11298422 :   x += 2;
    1722    11298422 :   z += 2;
    1723    11298422 :   y += 2;
    1724    11298422 :   gel(z,dz) = gcopy(p2);
    1725             : 
    1726    27765220 :   for (i=dx-1; i>=dy; i--)
    1727             :   {
    1728    16466798 :     av1=avma; p1=gel(x,i);
    1729    16466798 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1730    16466798 :     if (y_lead) p1 = simplify(f(p1,y_lead));
    1731             : 
    1732    16466798 :     if (isrationalzero(p1)) { set_avma(av1); p1 = gen_0; }
    1733             :     else
    1734     6921881 :       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
    1735    16466798 :     gel(z,i-dy) = p1;
    1736             :   }
    1737    11298422 :   if (!pr) return gerepileupto(av,z-2);
    1738             : 
    1739     4732902 :   rem = (GEN)avma; av1 = (pari_sp)new_chunk(dx+3);
    1740     5078723 :   for (sx=0; ; i--)
    1741             :   {
    1742     5424544 :     p1 = gel(x,i);
    1743             :     /* we always enter this loop at least once */
    1744     5078723 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1745     5078723 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1746     5078723 :     if (!gequal0(p1)) { sx = 1; break; } /* remainder is non-zero */
    1747     3414092 :     if (!isexactzero(p1)) break;
    1748     3401157 :     if (!i) break;
    1749      345821 :     set_avma(av1);
    1750             :   }
    1751     4732902 :   if (pr == ONLY_DIVIDES)
    1752             :   {
    1753        1568 :     if (sx) return gc_NULL(av);
    1754        1561 :     avma = (pari_sp)rem;
    1755        1561 :     return gerepileupto(av,z-2);
    1756             :   }
    1757     4731334 :   lr=i+3; rem -= lr;
    1758     4731334 :   if (avma==av1) { avma = (pari_sp)rem; p1 = gcopy(p1); }
    1759     4688198 :   else p1 = gerepileupto((pari_sp)rem,p1);
    1760     4731334 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1761     4731334 :   rem[1] = z[-1];
    1762     4731334 :   rem += 2;
    1763     4731334 :   gel(rem,i) = p1;
    1764     5881674 :   for (i--; i>=0; i--)
    1765             :   {
    1766     1150340 :     av1=avma; p1 = gel(x,i);
    1767     1150340 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1768     1150340 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1769     1150340 :     gel(rem,i) = avma==av1? gcopy(p1):gerepileupto(av1,p1);
    1770             :   }
    1771     4731334 :   rem -= 2;
    1772     4731334 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1773     4731334 :   if (pr == ONLY_REM) return gerepileupto(av,rem);
    1774     4071487 :   z -= 2;
    1775             :   {
    1776     4071487 :     GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
    1777     4071487 :     gerepilemanysp(av,avy,gptr,2); *pr = rem; return z;
    1778             :   }
    1779             : }
    1780             : 
    1781             : GEN
    1782    11288835 : RgX_divrem(GEN x, GEN y, GEN *pr)
    1783             : {
    1784    11288835 :   if (pr == ONLY_REM) return RgX_rem(x, y);
    1785    11288835 :   return RgX_divrem_i(x, y, pr);
    1786             : }
    1787             : 
    1788             : /* x and y in (R[Y]/T)[X]  (lifted), T in R[Y]. y preferably monic */
    1789             : GEN
    1790       69350 : RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)
    1791             : {
    1792             :   long vx, dx, dy, dz, i, j, sx, lr;
    1793             :   pari_sp av0, av, tetpil;
    1794             :   GEN z,p1,rem,lead;
    1795             : 
    1796       69350 :   if (!signe(y)) pari_err_INV("RgXQX_divrem",y);
    1797       69350 :   vx = varn(x);
    1798       69350 :   dx = degpol(x);
    1799       69350 :   dy = degpol(y);
    1800       69350 :   if (dx < dy)
    1801             :   {
    1802       22869 :     if (pr)
    1803             :     {
    1804       22869 :       av0 = avma; x = RgXQX_red(x, T);
    1805       22869 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: gen_0; }
    1806       22869 :       if (pr == ONLY_REM) return x;
    1807           0 :       *pr = x;
    1808             :     }
    1809           0 :     return pol_0(vx);
    1810             :   }
    1811       46481 :   lead = leading_coeff(y);
    1812       46481 :   if (!dy) /* y is constant */
    1813             :   {
    1814         546 :     if (pr && pr != ONLY_DIVIDES)
    1815             :     {
    1816           0 :       if (pr == ONLY_REM) return pol_0(vx);
    1817           0 :       *pr = pol_0(vx);
    1818             :     }
    1819         546 :     if (gequal1(lead)) return RgX_copy(x);
    1820           0 :     av0 = avma; x = gmul(x, ginvmod(lead,T)); tetpil = avma;
    1821           0 :     return gerepile(av0,tetpil,RgXQX_red(x,T));
    1822             :   }
    1823       45935 :   av0 = avma; dz = dx-dy;
    1824       45935 :   lead = gequal1(lead)? NULL: gclone(ginvmod(lead,T));
    1825       45935 :   set_avma(av0);
    1826       45935 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1827       45935 :   x += 2; y += 2; z += 2;
    1828             : 
    1829       45935 :   p1 = gel(x,dx); av = avma;
    1830       45935 :   gel(z,dz) = lead? gerepileupto(av, grem(gmul(p1,lead), T)): gcopy(p1);
    1831      288820 :   for (i=dx-1; i>=dy; i--)
    1832             :   {
    1833      242885 :     av=avma; p1=gel(x,i);
    1834      242885 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1835      242885 :     if (lead) p1 = gmul(grem(p1, T), lead);
    1836      242885 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil, grem(p1, T));
    1837             :   }
    1838       45935 :   if (!pr) { if (lead) gunclone(lead); return z-2; }
    1839             : 
    1840       44689 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    1841       70966 :   for (sx=0; ; i--)
    1842             :   {
    1843       97243 :     p1 = gel(x,i);
    1844       70966 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1845       70966 :     tetpil=avma; p1 = grem(p1, T); if (!gequal0(p1)) { sx = 1; break; }
    1846       31032 :     if (!i) break;
    1847       26277 :     set_avma(av);
    1848             :   }
    1849       44689 :   if (pr == ONLY_DIVIDES)
    1850             :   {
    1851        2305 :     if (lead) gunclone(lead);
    1852        2305 :     if (sx) return gc_NULL(av0);
    1853        2207 :     avma = (pari_sp)rem; return z-2;
    1854             :   }
    1855       42384 :   lr=i+3; rem -= lr;
    1856       42384 :   rem[0] = evaltyp(t_POL) | evallg(lr);
    1857       42384 :   rem[1] = z[-1];
    1858       42384 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    1859       42384 :   rem += 2; gel(rem,i) = p1;
    1860      115981 :   for (i--; i>=0; i--)
    1861             :   {
    1862       73597 :     av=avma; p1 = gel(x,i);
    1863      273148 :     for (j=0; j<=i && j<=dz; j++)
    1864      199551 :       p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1865       73597 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, grem(p1, T));
    1866             :   }
    1867       42384 :   rem -= 2;
    1868       42384 :   if (lead) gunclone(lead);
    1869       42384 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1870       42384 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    1871          91 :   *pr = rem; return z-2;
    1872             : }
    1873             : 
    1874             : /*******************************************************************/
    1875             : /*                                                                 */
    1876             : /*                        PSEUDO-DIVISION                          */
    1877             : /*                                                                 */
    1878             : /*******************************************************************/
    1879             : INLINE GEN
    1880      462189 : rem(GEN c, GEN T)
    1881             : {
    1882      462189 :   if (T && typ(c) == t_POL && varn(c) == varn(T)) c = RgX_rem(c, T);
    1883      462189 :   return c;
    1884             : }
    1885             : 
    1886             : /* x, y, are ZYX, lc(y) is an integer, T is a ZY */
    1887             : int
    1888        6136 : ZXQX_dvd(GEN x, GEN y, GEN T)
    1889             : {
    1890             :   long dx, dy, dz, i, p, T_ismonic;
    1891        6136 :   pari_sp av = avma, av2;
    1892             :   GEN y_lead;
    1893             : 
    1894        6136 :   if (!signe(y)) pari_err_INV("ZXQX_dvd",y);
    1895        6136 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1896        6136 :   if (typ(y_lead) == t_POL) y_lead = gel(y_lead, 2); /* t_INT */
    1897             :   /* if monic, no point in using pseudo-division */
    1898        6136 :   if (gequal1(y_lead)) return signe(RgXQX_rem(x, y, T)) == 0;
    1899        3735 :   T_ismonic = gequal1(leading_coeff(T));
    1900        3735 :   dx = degpol(x);
    1901        3735 :   if (dx < dy) return !signe(x);
    1902        3735 :   (void)new_chunk(2);
    1903        3735 :   x = RgX_recip_shallow(x)+2;
    1904        3735 :   y = RgX_recip_shallow(y)+2;
    1905             :   /* pay attention to sparse divisors */
    1906        7652 :   for (i = 1; i <= dy; i++)
    1907        3917 :     if (!signe(gel(y,i))) gel(y,i) = NULL;
    1908        3735 :   dz = dx-dy; p = dz+1;
    1909        3735 :   av2 = avma;
    1910             :   for (;;)
    1911       30634 :   {
    1912       34369 :     GEN m, x0 = gel(x,0), y0 = y_lead, cx = content(x0);
    1913       34369 :     x0 = gneg(x0); p--;
    1914       34369 :     m = gcdii(cx, y0);
    1915       34369 :     if (!equali1(m))
    1916             :     {
    1917       33319 :       x0 = gdiv(x0, m);
    1918       33319 :       y0 = diviiexact(y0, m);
    1919       33319 :       if (equali1(y0)) y0 = NULL;
    1920             :     }
    1921       69928 :     for (i=1; i<=dy; i++)
    1922             :     {
    1923       35559 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1924       35559 :       if (gel(y,i)) c = gadd(c, gmul(x0,gel(y,i)));
    1925       35559 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1926       35559 :       gel(x,i) = c;
    1927             :     }
    1928      399662 :     for (   ; i<=dx; i++)
    1929             :     {
    1930      365293 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    1931      365293 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    1932      365293 :       gel(x,i) = c;
    1933             :     }
    1934       38097 :     do { x++; dx--; } while (dx >= 0 && !signe(gel(x,0)));
    1935       34369 :     if (dx < dy) break;
    1936       30634 :     if (gc_needed(av2,1))
    1937             :     {
    1938           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZXQX_dvd dx = %ld >= %ld",dx,dy);
    1939           0 :       gerepilecoeffs(av2,x,dx+1);
    1940             :     }
    1941             :   }
    1942        3735 :   return gc_bool(av, dx < 0);
    1943             : }
    1944             : 
    1945             : /* T either NULL or a t_POL. */
    1946             : GEN
    1947       27759 : RgXQX_pseudorem(GEN x, GEN y, GEN T)
    1948             : {
    1949       27759 :   long vx = varn(x), dx, dy, dz, i, lx, p;
    1950       27759 :   pari_sp av = avma, av2;
    1951             :   GEN y_lead;
    1952             : 
    1953       27759 :   if (!signe(y)) pari_err_INV("RgXQX_pseudorem",y);
    1954       27759 :   dy = degpol(y); y_lead = gel(y,dy+2);
    1955             :   /* if monic, no point in using pseudo-division */
    1956       27759 :   if (gequal1(y_lead)) return T? RgXQX_rem(x, y, T): RgX_rem(x, y);
    1957       24490 :   dx = degpol(x);
    1958       24490 :   if (dx < dy) return RgX_copy(x);
    1959       24490 :   (void)new_chunk(2);
    1960       24490 :   x = RgX_recip_shallow(x)+2;
    1961       24490 :   y = RgX_recip_shallow(y)+2;
    1962             :   /* pay attention to sparse divisors */
    1963       74723 :   for (i = 1; i <= dy; i++)
    1964       50233 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    1965       24490 :   dz = dx-dy; p = dz+1;
    1966       24490 :   av2 = avma;
    1967             :   for (;;)
    1968             :   {
    1969      148286 :     gel(x,0) = gneg(gel(x,0)); p--;
    1970      273384 :     for (i=1; i<=dy; i++)
    1971             :     {
    1972      186996 :       GEN c = gmul(y_lead, gel(x,i));
    1973      186996 :       if (gel(y,i)) c = gadd(c, gmul(gel(x,0),gel(y,i)));
    1974      186996 :       gel(x,i) = rem(c, T);
    1975             :     }
    1976      295415 :     for (   ; i<=dx; i++)
    1977             :     {
    1978      209027 :       GEN c = gmul(y_lead, gel(x,i));
    1979      209027 :       gel(x,i) = rem(c, T);
    1980             :     }
    1981       95509 :     do { x++; dx--; } while (dx >= 0 && gequal0(gel(x,0)));
    1982       86388 :     if (dx < dy) break;
    1983       61898 :     if (gc_needed(av2,1))
    1984             :     {
    1985           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudorem dx = %ld >= %ld",dx,dy);
    1986           0 :       gerepilecoeffs(av2,x,dx+1);
    1987             :     }
    1988             :   }
    1989       24490 :   if (dx < 0) return pol_0(vx);
    1990       22453 :   lx = dx+3; x -= 2;
    1991       22453 :   x[0] = evaltyp(t_POL) | evallg(lx);
    1992       22453 :   x[1] = evalsigne(1) | evalvarn(vx);
    1993       22453 :   x = RgX_recip_shallow(x);
    1994       22453 :   if (p)
    1995             :   { /* multiply by y[0]^p   [beware dummy vars from FpX_FpXY_resultant] */
    1996        2506 :     GEN t = y_lead;
    1997        2506 :     if (T && typ(t) == t_POL && varn(t) == varn(T))
    1998           0 :       t = RgXQ_powu(t, p, T);
    1999             :     else
    2000        2506 :       t = gpowgs(t, p);
    2001        6825 :     for (i=2; i<lx; i++)
    2002             :     {
    2003        4319 :       GEN c = gmul(gel(x,i), t);
    2004        4319 :       gel(x,i) = rem(c,T);
    2005             :     }
    2006        2506 :     if (!T) return gerepileupto(av, x);
    2007             :   }
    2008       19947 :   return gerepilecopy(av, x);
    2009             : }
    2010             : 
    2011             : GEN
    2012       27759 : RgX_pseudorem(GEN x, GEN y) { return RgXQX_pseudorem(x,y, NULL); }
    2013             : 
    2014             : /* Compute z,r s.t lc(y)^(dx-dy+1) x = z y + r */
    2015             : GEN
    2016        5333 : RgXQX_pseudodivrem(GEN x, GEN y, GEN T, GEN *ptr)
    2017             : {
    2018        5333 :   long vx = varn(x), dx, dy, dz, i, iz, lx, lz, p;
    2019        5333 :   pari_sp av = avma, av2;
    2020             :   GEN z, r, ypow, y_lead;
    2021             : 
    2022        5333 :   if (!signe(y)) pari_err_INV("RgXQX_pseudodivrem",y);
    2023        5333 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2024        5333 :   if (gequal1(y_lead)) return T? RgXQX_divrem(x,y, T, ptr): RgX_divrem(x,y, ptr);
    2025        4780 :   dx = degpol(x);
    2026        4780 :   if (dx < dy) { *ptr = RgX_copy(x); return pol_0(vx); }
    2027        4780 :   if (dx == dy)
    2028             :   {
    2029          28 :     GEN x_lead = gel(x,lg(x)-1);
    2030          28 :     x = RgX_renormalize_lg(leafcopy(x), lg(x)-1);
    2031          28 :     y = RgX_renormalize_lg(leafcopy(y), lg(y)-1);
    2032          28 :     r = RgX_sub(RgX_Rg_mul(x, y_lead), RgX_Rg_mul(y, x_lead));
    2033          28 :     *ptr = gerepileupto(av, r); return scalarpol(x_lead, vx);
    2034             :   }
    2035        4752 :   (void)new_chunk(2);
    2036        4752 :   x = RgX_recip_shallow(x)+2;
    2037        4752 :   y = RgX_recip_shallow(y)+2;
    2038             :   /* pay attention to sparse divisors */
    2039       23489 :   for (i = 1; i <= dy; i++)
    2040       18737 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2041        4752 :   dz = dx-dy; p = dz+1;
    2042        4752 :   lz = dz+3;
    2043        4752 :   z = cgetg(lz, t_POL);
    2044        4752 :   z[1] = evalsigne(1) | evalvarn(vx);
    2045        4752 :   for (i = 2; i < lz; i++) gel(z,i) = gen_0;
    2046        4752 :   ypow = new_chunk(dz+1);
    2047        4752 :   gel(ypow,0) = gen_1;
    2048        4752 :   gel(ypow,1) = y_lead;
    2049        8028 :   for (i=2; i<=dz; i++)
    2050             :   {
    2051        3276 :     GEN c = gmul(gel(ypow,i-1), y_lead);
    2052        3276 :     gel(ypow,i) = rem(c,T);
    2053             :   }
    2054        4752 :   av2 = avma;
    2055        4752 :   for (iz=2;;)
    2056             :   {
    2057       15400 :     p--;
    2058       10076 :     gel(z,iz++) = rem(gmul(gel(x,0), gel(ypow,p)), T);
    2059       44531 :     for (i=1; i<=dy; i++)
    2060             :     {
    2061       34455 :       GEN c = gmul(y_lead, gel(x,i));
    2062       34455 :       if (gel(y,i)) c = gsub(c, gmul(gel(x,0),gel(y,i)));
    2063       34455 :       gel(x,i) = rem(c, T);
    2064             :     }
    2065       24116 :     for (   ; i<=dx; i++)
    2066             :     {
    2067       14040 :       GEN c = gmul(y_lead, gel(x,i));
    2068       14040 :       gel(x,i) = rem(c,T);
    2069             :     }
    2070       10076 :     x++; dx--;
    2071       10076 :     while (dx >= dy && gequal0(gel(x,0))) { x++; dx--; iz++; }
    2072       10076 :     if (dx < dy) break;
    2073        5324 :     if (gc_needed(av2,1))
    2074             :     {
    2075           0 :       GEN X = x-2;
    2076           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudodivrem dx=%ld >= %ld",dx,dy);
    2077           0 :       X[0] = evaltyp(t_POL)|evallg(dx+3); X[1] = z[1]; /* hack */
    2078           0 :       gerepileall(av2,2, &X, &z); x = X+2;
    2079             :     }
    2080             :   }
    2081        4752 :   while (dx >= 0 && gequal0(gel(x,0))) { x++; dx--; }
    2082        4752 :   if (dx < 0)
    2083          98 :     x = pol_0(vx);
    2084             :   else
    2085             :   {
    2086        4654 :     lx = dx+3; x -= 2;
    2087        4654 :     x[0] = evaltyp(t_POL) | evallg(lx);
    2088        4654 :     x[1] = evalsigne(1) | evalvarn(vx);
    2089        4654 :     x = RgX_recip_shallow(x);
    2090             :   }
    2091        4752 :   z = RgX_recip_shallow(z);
    2092        4752 :   r = x;
    2093        4752 :   if (p)
    2094             :   {
    2095        1674 :     GEN c = gel(ypow,p); r = RgX_Rg_mul(r, c);
    2096        1674 :     if (T && typ(c) == t_POL && varn(c) == varn(T)) r = RgXQX_red(r, T);
    2097             :   }
    2098        4752 :   gerepileall(av, 2, &z, &r);
    2099        4752 :   *ptr = r; return z;
    2100             : }
    2101             : GEN
    2102        5165 : RgX_pseudodivrem(GEN x, GEN y, GEN *ptr)
    2103        5165 : { return RgXQX_pseudodivrem(x,y,NULL,ptr); }
    2104             : 
    2105             : GEN
    2106       12789 : RgXQX_mul(GEN x, GEN y, GEN T)
    2107             : {
    2108       12789 :   return RgXQX_red(RgX_mul(x,y), T);
    2109             : }
    2110             : GEN
    2111    63340243 : RgX_Rg_mul(GEN y, GEN x) {
    2112             :   long i, ly;
    2113    63340243 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2114    63340243 :   if (ly == 2) return z;
    2115    63255515 :   for (i = 2; i < ly; i++) gel(z,i) = gmul(x,gel(y,i));
    2116    63255508 :   return normalizepol_lg(z,ly);
    2117             : }
    2118             : GEN
    2119       14427 : RgX_muls(GEN y, long x) {
    2120             :   long i, ly;
    2121       14427 :   GEN z = cgetg_copy(y, &ly); z[1] = y[1];
    2122       14427 :   if (ly == 2) return z;
    2123       14392 :   for (i = 2; i < ly; i++) gel(z,i) = gmulsg(x,gel(y,i));
    2124       14392 :   return normalizepol_lg(z,ly);
    2125             : }
    2126             : GEN
    2127        4613 : RgXQX_RgXQ_mul(GEN x, GEN y, GEN T)
    2128             : {
    2129        4613 :   return RgXQX_red(RgX_Rg_mul(x,y), T);
    2130             : }
    2131             : GEN
    2132          77 : RgXQV_RgXQ_mul(GEN v, GEN x, GEN T)
    2133             : {
    2134          77 :   return RgXQV_red(RgV_Rg_mul(v,x), T);
    2135             : }
    2136             : 
    2137             : GEN
    2138        1792 : RgXQX_sqr(GEN x, GEN T)
    2139             : {
    2140        1792 :   return RgXQX_red(RgX_sqr(x), T);
    2141             : }
    2142             : 
    2143             : GEN
    2144         448 : RgXQX_powers(GEN P, long n, GEN T)
    2145             : {
    2146         448 :   GEN v = cgetg(n+2, t_VEC);
    2147             :   long i;
    2148         448 :   gel(v, 1) = pol_1(varn(T));
    2149         448 :   if (n==0) return v;
    2150         448 :   gel(v, 2) = gcopy(P);
    2151         448 :   for (i = 2; i <= n; i++) gel(v,i+1) = RgXQX_mul(P, gel(v,i), T);
    2152         448 :   return v;
    2153             : }
    2154             : 
    2155             : static GEN
    2156      124397 : _add(void *data, GEN x, GEN y) { (void)data; return RgX_add(x, y); }
    2157             : static GEN
    2158           0 : _sub(void *data, GEN x, GEN y) { (void)data; return RgX_sub(x, y); }
    2159             : static GEN
    2160       55817 : _sqr(void *data, GEN x) { return RgXQ_sqr(x, (GEN)data); }
    2161             : static GEN
    2162       39663 : _mul(void *data, GEN x, GEN y) { return RgXQ_mul(x,y, (GEN)data); }
    2163             : static GEN
    2164      226555 : _cmul(void *data, GEN P, long a, GEN x) { (void)data; return RgX_Rg_mul(x,gel(P,a+2)); }
    2165             : static GEN
    2166      221676 : _one(void *data) { return pol_1(varn((GEN)data)); }
    2167             : static GEN
    2168         245 : _zero(void *data) { return pol_0(varn((GEN)data)); }
    2169             : static GEN
    2170      129948 : _red(void *data, GEN x) { (void)data; return gcopy(x); }
    2171             : 
    2172             : static struct bb_algebra RgXQ_algebra = { _red, _add, _sub,
    2173             :               _mul, _sqr, _one, _zero };
    2174             : 
    2175             : GEN
    2176           0 : RgX_RgXQV_eval(GEN Q, GEN x, GEN T)
    2177             : {
    2178           0 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)T,&RgXQ_algebra,_cmul);
    2179             : }
    2180             : 
    2181             : GEN
    2182      101682 : RgX_RgXQ_eval(GEN Q, GEN x, GEN T)
    2183             : {
    2184      101682 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2185      101682 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)T,&RgXQ_algebra,_cmul);
    2186             : }
    2187             : 
    2188             : /* mod X^n */
    2189             : struct modXn {
    2190             :   long v; /* varn(X) */
    2191             :   long n;
    2192             : } ;
    2193             : static GEN
    2194        2254 : _sqrXn(void *data, GEN x) {
    2195        2254 :   struct modXn *S = (struct modXn*)data;
    2196        2254 :   return RgXn_sqr(x, S->n);
    2197             : }
    2198             : static GEN
    2199        1295 : _mulXn(void *data, GEN x, GEN y) {
    2200        1295 :   struct modXn *S = (struct modXn*)data;
    2201        1295 :   return RgXn_mul(x,y, S->n);
    2202             : }
    2203             : static GEN
    2204        1477 : _oneXn(void *data) {
    2205        1477 :   struct modXn *S = (struct modXn*)data;
    2206        1477 :   return pol_1(S->v);
    2207             : }
    2208             : static GEN
    2209           0 : _zeroXn(void *data) {
    2210           0 :   struct modXn *S = (struct modXn*)data;
    2211           0 :   return pol_0(S->v);
    2212             : }
    2213             : static struct bb_algebra RgXn_algebra = { _red, _add, _sub, _mulXn, _sqrXn,
    2214             :                                           _oneXn, _zeroXn };
    2215             : 
    2216             : GEN
    2217         357 : RgXn_powers(GEN x, long m, long n)
    2218             : {
    2219         357 :   long d = degpol(x);
    2220         357 :   int use_sqr = (d<<1) >= n;
    2221             :   struct modXn S;
    2222         357 :   S.v = varn(x); S.n = n;
    2223         357 :   return gen_powers(x,m,use_sqr,(void*)&S,_sqrXn,_mulXn,_oneXn);
    2224             : }
    2225             : 
    2226             : GEN
    2227        1834 : RgXn_powu_i(GEN x, ulong m, long n)
    2228             : {
    2229             :   struct modXn S;
    2230        1834 :   S.v = varn(x); S.n = n;
    2231        1834 :   return gen_powu_i(x, m, (void*)&S,_sqrXn,_mulXn);
    2232             : }
    2233             : GEN
    2234           0 : RgXn_powu(GEN x, ulong m, long n)
    2235             : {
    2236             :   struct modXn S;
    2237           0 :   S.v = varn(x); S.n = n;
    2238           0 :   return gen_powu(x, m, (void*)&S,_sqrXn,_mulXn);
    2239             : }
    2240             : 
    2241             : GEN
    2242         714 : RgX_RgXnV_eval(GEN Q, GEN x, long n)
    2243             : {
    2244             :   struct modXn S;
    2245         714 :   S.v = varn(gel(x,2)); S.n = n;
    2246         714 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)&S,&RgXn_algebra,_cmul);
    2247             : }
    2248             : 
    2249             : GEN
    2250           0 : RgX_RgXn_eval(GEN Q, GEN x, long n)
    2251             : {
    2252           0 :   int use_sqr = 2*degpol(x) >= n;
    2253             :   struct modXn S;
    2254           0 :   S.v = varn(x); S.n = n;
    2255           0 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2256             : }
    2257             : 
    2258             : /* Q(x) mod t^n, x in R[t], n >= 1 */
    2259             : GEN
    2260        6503 : RgXn_eval(GEN Q, GEN x, long n)
    2261             : {
    2262        6503 :   long d = degpol(x);
    2263             :   int use_sqr;
    2264             :   struct modXn S;
    2265        6503 :   if (d == 1 && isrationalzero(gel(x,2)))
    2266             :   {
    2267        6496 :     GEN y = RgX_unscale(Q, gel(x,3));
    2268        6496 :     setvarn(y, varn(x)); return y;
    2269             :   }
    2270           7 :   S.v = varn(x);
    2271           7 :   S.n = n;
    2272           7 :   use_sqr = (d<<1) >= n;
    2273           7 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2274             : }
    2275             : 
    2276             : /* (f*g mod t^n) \ t^n2, assuming 2*n2 >= n */
    2277             : static GEN
    2278       43588 : RgXn_mulhigh(GEN f, GEN g, long n2, long n)
    2279             : {
    2280       43588 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2281       43588 :   return RgX_add(RgX_mulhigh_i(fl, g, n2), RgXn_mul(fh, g, n - n2));
    2282             : }
    2283             : 
    2284             : /* (f^2 mod t^n) \ t^n2, assuming 2*n2 >= n */
    2285             : static GEN
    2286           0 : RgXn_sqrhigh(GEN f, long n2, long n)
    2287             : {
    2288           0 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2289           0 :   return RgX_add(RgX_mulhigh_i(fl, f, n2), RgXn_mul(fh, f, n - n2));
    2290             : }
    2291             : 
    2292             : GEN
    2293      126154 : RgXn_inv_i(GEN f, long e)
    2294             : {
    2295             :   pari_sp av;
    2296             :   ulong mask;
    2297             :   GEN W, a;
    2298      126154 :   long v = varn(f), n = 1;
    2299             : 
    2300      126154 :   if (!signe(f)) pari_err_INV("RgXn_inv",f);
    2301      126154 :   a = ginv(gel(f,2));
    2302      126154 :   if (e == 1) return scalarpol(a, v);
    2303      126126 :   else if (e == 2)
    2304             :   {
    2305             :     GEN b;
    2306      108745 :     if (degpol(f) <= 0 || gequal0(b = gel(f,3))) return scalarpol(a, v);
    2307       88690 :     av = avma; b = gneg(b);
    2308       88690 :     if (!gequal1(a)) b = gmul(b, gsqr(a));
    2309       88690 :     W = deg1pol_shallow(b, a, v);
    2310       88690 :     return gcopy(W);
    2311             :   }
    2312       17381 :   W = scalarpol_shallow(a,v);
    2313       17381 :   mask = quadratic_prec_mask(e);
    2314       17381 :   av = avma;
    2315       78350 :   for (;mask>1;)
    2316             :   {
    2317             :     GEN u, fr;
    2318       43588 :     long n2 = n;
    2319       43588 :     n<<=1; if (mask & 1) n--;
    2320       43588 :     mask >>= 1;
    2321       43588 :     fr = RgXn_red_shallow(f, n);
    2322       43588 :     u = RgXn_mul(W, RgXn_mulhigh(fr, W, n2, n), n-n2);
    2323       43588 :     W = RgX_sub(W, RgX_shift_shallow(u, n2));
    2324       43588 :     if (gc_needed(av,2))
    2325             :     {
    2326           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_inv, e = %ld", n);
    2327           0 :       W = gerepileupto(av, W);
    2328             :     }
    2329             :   }
    2330       17381 :   return W;
    2331             : }
    2332             : 
    2333             : static GEN
    2334           0 : RgXn_inv_FpX(GEN x, long e, GEN p)
    2335             : {
    2336           0 :   pari_sp av = avma;
    2337             :   GEN r;
    2338           0 :   if (lgefint(p) == 3)
    2339             :   {
    2340           0 :     ulong pp = uel(p, 2);
    2341           0 :     r = Flx_to_ZX_inplace(Flxn_inv(RgX_to_Flx(x, pp), e, pp));
    2342             :   }
    2343             :   else
    2344           0 :     r = FpXn_inv(RgX_to_FpX(x, p), e, p);
    2345           0 :   return gerepileupto(av, FpX_to_mod(r, p));
    2346             : }
    2347             : 
    2348             : static GEN
    2349           0 : RgXn_inv_FpXQX(GEN x, long n, GEN pol, GEN p)
    2350             : {
    2351           0 :   pari_sp av = avma;
    2352           0 :   GEN r, T = RgX_to_FpX(pol, p);
    2353           0 :   if (signe(T) == 0) pari_err_OP("/", gen_1, x);
    2354           0 :   r = FpXQXn_inv(RgX_to_FpXQX(x, T, p), n, T, p);
    2355           0 :   return gerepileupto(av, FpXQX_to_mod(r, pol, p));
    2356             : }
    2357             : 
    2358             : #define code(t1,t2) ((t1 << 6) | t2)
    2359             : 
    2360             : static GEN
    2361       82957 : RgXn_inv_fast(GEN x, long e)
    2362             : {
    2363             :   GEN p, pol;
    2364             :   long pa;
    2365       82957 :   long t = RgX_type(x,&p,&pol,&pa);
    2366       82957 :   switch(t)
    2367             :   {
    2368           0 :     case t_INTMOD: return RgXn_inv_FpX(x, e, p);
    2369             :     case code(t_POLMOD, t_INTMOD):
    2370           0 :                    return RgXn_inv_FpXQX(x, e, pol, p);
    2371       82957 :     default:       return NULL;
    2372             :   }
    2373             : }
    2374             : #undef code
    2375             : 
    2376             : GEN
    2377       82957 : RgXn_inv(GEN f, long e)
    2378             : {
    2379       82957 :   GEN h = RgXn_inv_fast(f, e);
    2380       82957 :   if (h) return h;
    2381       82957 :   return RgXn_inv_i(f, e);
    2382             : }
    2383             : 
    2384             : GEN
    2385       12957 : RgXn_exp(GEN h, long e)
    2386             : {
    2387       12957 :   pari_sp av = avma, av2;
    2388       12957 :   long v = varn(h), n=1;
    2389       12957 :   GEN f = pol_1(v), g = pol_1(v);
    2390       12957 :   ulong mask = quadratic_prec_mask(e);
    2391       12957 :   av2 = avma;
    2392       12957 :   if (signe(h)==0 || degpol(h)<1 || !gequal0(gel(h,2)))
    2393           0 :     pari_err_DOMAIN("RgXn_exp","valuation", "<", gen_1, h);
    2394       40110 :   for (;mask>1;)
    2395             :   {
    2396             :     GEN q, w;
    2397       14196 :     long n2 = n;
    2398       14196 :     n<<=1; if (mask & 1) n--;
    2399       14196 :     mask >>= 1;
    2400       14196 :     g = RgX_sub(RgX_muls(g,2),RgXn_mul(f,RgXn_sqr(g,n2),n2));
    2401       14196 :     q = RgX_deriv(RgXn_red_shallow(h,n2));
    2402       14196 :     w = RgX_add(q, RgXn_mul(g, RgX_sub(RgX_deriv(f), RgXn_mul(f,q,n-1)),n-1));
    2403       14196 :     f = RgX_add(f, RgXn_mul(f, RgX_sub(RgXn_red_shallow(h, n), RgX_integ(w)), n));
    2404       14196 :     if (gc_needed(av2,2))
    2405             :     {
    2406           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_exp, e = %ld", n);
    2407           0 :       gerepileall(av2, 2, &f, &g);
    2408             :     }
    2409             :   }
    2410       12957 :   return gerepileupto(av, f);
    2411             : }
    2412             : 
    2413             : GEN
    2414          91 : RgXn_reverse(GEN f, long e)
    2415             : {
    2416          91 :   pari_sp av = avma, av2;
    2417             :   ulong mask;
    2418             :   GEN fi, a, df, W, an;
    2419          91 :   long v = varn(f), n=1;
    2420          91 :   if (degpol(f)<1 || !gequal0(gel(f,2)))
    2421           0 :     pari_err_INV("serreverse",f);
    2422          91 :   fi = ginv(gel(f,3));
    2423          91 :   a = deg1pol_shallow(fi,gen_0,v);
    2424          91 :   if (e <= 2) return gerepilecopy(av, a);
    2425          91 :   W = scalarpol(fi,v);
    2426          91 :   df = RgX_deriv(f);
    2427          91 :   mask = quadratic_prec_mask(e);
    2428          91 :   av2 = avma;
    2429         539 :   for (;mask>1;)
    2430             :   {
    2431             :     GEN u, fa, fr;
    2432         357 :     long n2 = n, rt;
    2433         357 :     n<<=1; if (mask & 1) n--;
    2434         357 :     mask >>= 1;
    2435         357 :     fr = RgXn_red_shallow(f, n);
    2436         357 :     rt = brent_kung_optpow(degpol(fr), 4, 3);
    2437         357 :     an = RgXn_powers(a, rt, n);
    2438         357 :     if (n>1)
    2439             :     {
    2440         357 :       long n4 = (n2+1)>>1;
    2441         357 :       GEN dfr = RgXn_red_shallow(df, n2);
    2442         357 :       dfr = RgX_RgXnV_eval(dfr, RgXnV_red_shallow(an, n2), n2);
    2443         357 :       u = RgX_shift(RgX_Rg_sub(RgXn_mul(W, dfr, n2), gen_1), -n4);
    2444         357 :       W = RgX_sub(W, RgX_shift(RgXn_mul(u, W, n2-n4), n4));
    2445             :     }
    2446         357 :     fa = RgX_sub(RgX_RgXnV_eval(fr, an, n), pol_x(v));
    2447         357 :     fa = RgX_shift(fa, -n2);
    2448         357 :     a = RgX_sub(a, RgX_shift(RgXn_mul(W, fa, n-n2), n2));
    2449         357 :     if (gc_needed(av2,2))
    2450             :     {
    2451           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_reverse, e = %ld", n);
    2452           0 :       gerepileall(av2, 2, &a, &W);
    2453             :     }
    2454             :   }
    2455          91 :   return gerepileupto(av, a);
    2456             : }
    2457             : 
    2458             : GEN
    2459           0 : RgXn_sqrt(GEN h, long e)
    2460             : {
    2461           0 :   pari_sp av = avma, av2;
    2462           0 :   long v = varn(h), n = 1;
    2463           0 :   GEN f = scalarpol(gen_1, v), df = f;
    2464           0 :   ulong mask = quadratic_prec_mask(e);
    2465           0 :   if (degpol(h)<0 || !gequal1(gel(h,2)))
    2466           0 :     pari_err_SQRTN("RgXn_sqrt",h);
    2467           0 :   av2 = avma;
    2468             :   while(1)
    2469           0 :   {
    2470           0 :     long n2 = n, m;
    2471             :     GEN g;
    2472           0 :     n<<=1; if (mask & 1) n--;
    2473           0 :     mask >>= 1;
    2474           0 :     m = n-n2;
    2475           0 :     g = RgX_sub(RgXn_sqrhigh(f, n2, n), RgX_shift_shallow(RgXn_red_shallow(h, n),-n2));
    2476           0 :     f = RgX_sub(f, RgX_shift_shallow(RgXn_mul(gmul2n(df, -1), g, m), n2));
    2477           0 :     if (mask==1) return gerepileupto(av, f);
    2478           0 :     g = RgXn_mul(df, RgXn_mulhigh(df, f, n2, n), m);
    2479           0 :     df = RgX_sub(df, RgX_shift_shallow(g, n2));
    2480           0 :     if (gc_needed(av2,2))
    2481             :     {
    2482           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_sqrt, e = %ld", n);
    2483           0 :       gerepileall(av2, 2, &f, &df);
    2484             :     }
    2485             :   }
    2486             : }
    2487             : 
    2488             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2489             : GEN
    2490       18186 : RgXQ_powu(GEN x, ulong n, GEN T)
    2491             : {
    2492             :   pari_sp av;
    2493             :   GEN y;
    2494             : 
    2495       18186 :   if (!n) return pol_1(varn(x));
    2496       16576 :   if (n == 1) return RgX_copy(x);
    2497       14868 :   av = avma; y = gen_powu(x, n, (void*)T, &_sqr, &_mul);
    2498       14868 :   return gerepileupto(av, y);
    2499             : }
    2500             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2501             : GEN
    2502       19354 : RgXQ_pow(GEN x, GEN n, GEN T)
    2503             : {
    2504             :   pari_sp av;
    2505       19354 :   long s = signe(n);
    2506             :   GEN y;
    2507             : 
    2508       19354 :   if (!s) return pol_1(varn(x));
    2509       19354 :   if (is_pm1(n) == 1)
    2510           0 :     return (s < 0)? RgXQ_inv(x, T): RgX_copy(x);
    2511       19354 :   av = avma;
    2512       19354 :   if (s < 0) x = RgXQ_inv(x, T);
    2513       19354 :   y = gen_pow(x, n, (void*)T, &_sqr, &_mul);
    2514       19354 :   return gerepileupto(av, y);
    2515             : }
    2516             : static GEN
    2517      158953 : _ZXQsqr(void *data, GEN x) { return ZXQ_sqr(x, (GEN)data); }
    2518             : static GEN
    2519       50429 : _ZXQmul(void *data, GEN x, GEN y) { return ZXQ_mul(x,y, (GEN)data); }
    2520             : /* x,T in Z[X], n in N, compute lift(x^n mod T)) */
    2521             : GEN
    2522      163569 : ZXQ_powu(GEN x, ulong n, GEN T)
    2523             : {
    2524             :   pari_sp av;
    2525             :   GEN y;
    2526             : 
    2527      163569 :   if (!n) return pol_1(varn(x));
    2528      163569 :   if (n == 1) return ZX_copy(x);
    2529      107343 :   av = avma; y = gen_powu(x, n, (void*)T, &_ZXQsqr, &_ZXQmul);
    2530      107352 :   return gerepileupto(av, y);
    2531             : }
    2532             : 
    2533             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2534             : GEN
    2535        3409 : RgXQ_powers(GEN x, long l, GEN T)
    2536             : {
    2537        3409 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2538        3409 :   return gen_powers(x, l, use_sqr, (void *)T,_sqr,_mul,_one);
    2539             : }
    2540             : 
    2541             : /* a in K = Q[X]/(T), returns [a^0, ..., a^n] */
    2542             : GEN
    2543        2338 : QXQ_powers(GEN a, long n, GEN T)
    2544             : {
    2545        2338 :   GEN den, v = RgXQ_powers(Q_remove_denom(a, &den), n, T);
    2546             :   /* den*a integral; v[i+1] = (den*a)^i in K */
    2547        2338 :   if (den)
    2548             :   { /* restore denominators */
    2549        1393 :     GEN d = den;
    2550             :     long i;
    2551        1393 :     gel(v,2) = a;
    2552        4109 :     for (i=3; i<=n+1; i++) {
    2553        2716 :       d = mulii(d,den);
    2554        2716 :       gel(v,i) = RgX_Rg_div(gel(v,i), d);
    2555             :     }
    2556             :   }
    2557        2338 :   return v;
    2558             : }
    2559             : 
    2560             : static GEN
    2561        1596 : do_QXQ_eval(GEN v, long imin, GEN a, GEN T)
    2562             : {
    2563        1596 :   long l, i, m = 0;
    2564             :   GEN dz, z;
    2565        1596 :   GEN V = cgetg_copy(v, &l);
    2566        5089 :   for (i = imin; i < l; i++)
    2567             :   {
    2568        3493 :     GEN c = gel(v, i);
    2569        3493 :     if (typ(c) == t_POL) m = maxss(m, degpol(c));
    2570             :   }
    2571        1596 :   z = Q_remove_denom(QXQ_powers(a, m, T), &dz);
    2572        1596 :   for (i = 1; i < imin; i++) V[i] = v[i];
    2573        5089 :   for (i = imin; i < l; i++)
    2574             :   {
    2575        3493 :     GEN c = gel(v,i);
    2576        3493 :     if (typ(c) == t_POL) c = QX_ZXQV_eval(c, z, dz);
    2577        3493 :     gel(V,i) = c;
    2578             :   }
    2579        1596 :   return V;
    2580             : }
    2581             : /* [ s(a mod T) | s <- lift(v) ], a,T are QX, v a QXV */
    2582             : GEN
    2583        1533 : QXV_QXQ_eval(GEN v, GEN a, GEN T)
    2584        1533 : { return do_QXQ_eval(v, 1, a, T); }
    2585             : GEN
    2586          63 : QXX_QXQ_eval(GEN v, GEN a, GEN T)
    2587          63 : { return normalizepol(do_QXQ_eval(v, 2, a, T)); }
    2588             : 
    2589             : GEN
    2590         287 : RgXQ_matrix_pow(GEN y, long n, long m, GEN P)
    2591             : {
    2592         287 :   return RgXV_to_RgM(RgXQ_powers(y,m-1,P),n);
    2593             : }
    2594             : 
    2595             : GEN
    2596       38668 : RgXQ_norm(GEN x, GEN T)
    2597             : {
    2598             :   pari_sp av;
    2599       38668 :   long dx = degpol(x);
    2600             :   GEN L, y;
    2601             : 
    2602       38668 :   av = avma; y = resultant(T, x);
    2603       38668 :   L = leading_coeff(T);
    2604       38668 :   if (gequal1(L) || !signe(x)) return y;
    2605           0 :   return gerepileupto(av, gdiv(y, gpowgs(L, dx)));
    2606             : }
    2607             : 
    2608             : GEN
    2609      113700 : RgX_blocks(GEN P, long n, long m)
    2610             : {
    2611      113700 :   GEN z = cgetg(m+1,t_VEC);
    2612      113700 :   long i,j, k=2, l = lg(P);
    2613      534272 :   for(i=1; i<=m; i++)
    2614             :   {
    2615      420572 :     GEN zi = cgetg(n+2,t_POL);
    2616      420572 :     zi[1] = P[1];
    2617      420572 :     gel(z,i) = zi;
    2618     2663251 :     for(j=2; j<n+2; j++)
    2619     2242679 :       gel(zi, j) = k==l ? gen_0 : gel(P,k++);
    2620      420572 :     zi = RgX_renormalize_lg(zi, n+2);
    2621             :   }
    2622      113700 :   return z;
    2623             : }
    2624             : 
    2625             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
    2626             : void
    2627       76936 : RgX_even_odd(GEN p, GEN *pe, GEN *po)
    2628             : {
    2629       76936 :   long n = degpol(p), v = varn(p), n0, n1, i;
    2630             :   GEN p0, p1;
    2631             : 
    2632       76938 :   if (n <= 0) { *pe = RgX_copy(p); *po = zeropol(v); return; }
    2633             : 
    2634       76651 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
    2635       76651 :   p0 = cgetg(n0+2, t_POL); p0[1] = evalvarn(v)|evalsigne(1);
    2636       76652 :   p1 = cgetg(n1+2, t_POL); p1[1] = evalvarn(v)|evalsigne(1);
    2637     4828642 :   for (i=0; i<n1; i++)
    2638             :   {
    2639     4751987 :     p0[2+i] = p[2+(i<<1)];
    2640     4751987 :     p1[2+i] = p[3+(i<<1)];
    2641             :   }
    2642       76655 :   if (n1 != n0)
    2643       42381 :     p0[2+i] = p[2+(i<<1)];
    2644       76655 :   *pe = normalizepol(p0);
    2645       76652 :   *po = normalizepol(p1);
    2646             : }
    2647             : 
    2648             : /* write p(X) = a_0(X^k) + Xa_1(X^k) + ... + X^(k-1)a_{k-1}(X^k), shallow function */
    2649             : GEN
    2650       40670 : RgX_splitting(GEN p, long k)
    2651             : {
    2652       40670 :   long n = degpol(p), v = varn(p), m, i, j, l;
    2653             :   GEN r;
    2654             : 
    2655       40670 :   m = n/k;
    2656       40670 :   r = cgetg(k+1,t_VEC);
    2657      224154 :   for(i=1; i<=k; i++)
    2658             :   {
    2659      183484 :     gel(r,i) = cgetg(m+3, t_POL);
    2660      183484 :     mael(r,i,1) = evalvarn(v)|evalsigne(1);
    2661             :   }
    2662      553245 :   for (j=1, i=0, l=2; i<=n; i++)
    2663             :   {
    2664      512575 :     gmael(r,j,l) = gel(p,2+i);
    2665      512575 :     if (j==k) { j=1; l++; } else j++;
    2666             :   }
    2667      224154 :   for(i=1; i<=k; i++)
    2668      183484 :     gel(r,i) = normalizepol_lg(gel(r,i),i<j?l+1:l);
    2669       40670 :   return r;
    2670             : }
    2671             : 
    2672             : /*******************************************************************/
    2673             : /*                                                                 */
    2674             : /*                        Kronecker form                           */
    2675             : /*                                                                 */
    2676             : /*******************************************************************/
    2677             : 
    2678             : /* z in R[Y] representing an elt in R[X,Y] mod T(Y) in Kronecker form,
    2679             :  * i.e subst(lift(z), x, y^(2deg(z)-1)). Recover the "real" z, with
    2680             :  * normalized coefficients */
    2681             : GEN
    2682      132048 : Kronecker_to_mod(GEN z, GEN T)
    2683             : {
    2684      132048 :   long i,j,lx,l = lg(z), N = (degpol(T)<<1) + 1;
    2685      132048 :   GEN x, t = cgetg(N,t_POL);
    2686      132048 :   t[1] = T[1];
    2687      132048 :   lx = (l-2) / (N-2); x = cgetg(lx+3,t_POL);
    2688      132048 :   x[1] = z[1];
    2689      132048 :   T = RgX_copy(T);
    2690     1195242 :   for (i=2; i<lx+2; i++, z+= N-2)
    2691             :   {
    2692     1063194 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
    2693     1063194 :     gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2694             :   }
    2695      132048 :   N = (l-2) % (N-2) + 2;
    2696      132048 :   for (j=2; j<N; j++) t[j] = z[j];
    2697      132048 :   gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2698      132048 :   return normalizepol_lg(x, i+1);
    2699             : }
    2700             : 
    2701             : /*******************************************************************/
    2702             : /*                                                                 */
    2703             : /*                        Domain detection                         */
    2704             : /*                                                                 */
    2705             : /*******************************************************************/
    2706             : 
    2707             : static GEN
    2708      307532 : zero_FpX_mod(GEN p, long v)
    2709             : {
    2710      307532 :   GEN r = cgetg(3,t_POL);
    2711      307532 :   r[1] = evalvarn(v);
    2712      307532 :   gel(r,2) = mkintmod(gen_0, icopy(p));
    2713      307532 :   return r;
    2714             : }
    2715             : 
    2716             : static GEN
    2717      494605 : RgX_mul_FpX(GEN x, GEN y, GEN p)
    2718             : {
    2719      494605 :   pari_sp av = avma;
    2720             :   GEN r;
    2721      494605 :   if (lgefint(p) == 3)
    2722             :   {
    2723      466226 :     ulong pp = uel(p, 2);
    2724      466226 :     r = Flx_to_ZX_inplace(Flx_mul(RgX_to_Flx(x, pp),
    2725             :                                   RgX_to_Flx(y, pp), pp));
    2726             :   }
    2727             :   else
    2728       28379 :     r = FpX_mul(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2729      494605 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2730      384795 :   return gerepileupto(av, FpX_to_mod(r, p));
    2731             : }
    2732             : 
    2733             : static GEN
    2734          21 : zero_FpXQX_mod(GEN pol, GEN p, long v)
    2735             : {
    2736          21 :   GEN r = cgetg(3,t_POL);
    2737          21 :   r[1] = evalvarn(v);
    2738          21 :   gel(r,2) = mkpolmod(mkintmod(gen_0, icopy(p)), gcopy(pol));
    2739          21 :   return r;
    2740             : }
    2741             : 
    2742             : static GEN
    2743          77 : RgX_mul_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2744             : {
    2745          77 :   pari_sp av = avma;
    2746             :   long dT;
    2747             :   GEN kx, ky, r;
    2748          77 :   GEN T = RgX_to_FpX(pol, p);
    2749          77 :   if (signe(T)==0) pari_err_OP("*", x, y);
    2750          70 :   dT = degpol(T);
    2751          70 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2752          70 :   ky = ZXX_to_Kronecker(RgX_to_FpXQX(y, T, p), dT);
    2753          70 :   r = FpX_mul(kx, ky, p);
    2754          70 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2755          56 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2756             : }
    2757             : 
    2758             : static GEN
    2759      335184 : RgX_liftred(GEN x, GEN T)
    2760      335184 : { return RgXQX_red(liftpol_shallow(x), T); }
    2761             : 
    2762             : static GEN
    2763      128077 : RgX_mul_QXQX(GEN x, GEN y, GEN T)
    2764             : {
    2765      128077 :   pari_sp av = avma;
    2766      128077 :   long dT = degpol(T);
    2767      128077 :   GEN r = QX_mul(ZXX_to_Kronecker(RgX_liftred(x, T), dT),
    2768             :                  ZXX_to_Kronecker(RgX_liftred(y, T), dT));
    2769      128077 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2770             : }
    2771             : 
    2772             : static GEN
    2773        1260 : RgX_sqr_FpX(GEN x, GEN p)
    2774             : {
    2775        1260 :   pari_sp av = avma;
    2776             :   GEN r;
    2777        1260 :   if (lgefint(p) == 3)
    2778             :   {
    2779        1056 :     ulong pp = uel(p, 2);
    2780        1056 :     r = Flx_to_ZX_inplace(Flx_sqr(RgX_to_Flx(x, pp), pp));
    2781             :   }
    2782             :   else
    2783         204 :     r = FpX_sqr(RgX_to_FpX(x, p), p);
    2784        1260 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2785        1162 :   return gerepileupto(av, FpX_to_mod(r, p));
    2786             : }
    2787             : 
    2788             : static GEN
    2789         154 : RgX_sqr_FpXQX(GEN x, GEN pol, GEN p)
    2790             : {
    2791         154 :   pari_sp av = avma;
    2792             :   long dT;
    2793         154 :   GEN kx, r, T = RgX_to_FpX(pol, p);
    2794         154 :   if (signe(T)==0) pari_err_OP("*",x,x);
    2795         147 :   dT = degpol(T);
    2796         147 :   kx = ZXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    2797         147 :   r = FpX_sqr(kx, p);
    2798         147 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2799         147 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    2800             : }
    2801             : 
    2802             : static GEN
    2803        3768 : RgX_sqr_QXQX(GEN x, GEN T)
    2804             : {
    2805        3768 :   pari_sp av = avma;
    2806        3768 :   long dT = degpol(T);
    2807        3768 :   GEN r = QX_sqr(ZXX_to_Kronecker(RgX_liftred(x, T), dT));
    2808        3768 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    2809             : }
    2810             : 
    2811             : static GEN
    2812      649362 : RgX_rem_FpX(GEN x, GEN y, GEN p)
    2813             : {
    2814      649362 :   pari_sp av = avma;
    2815             :   GEN r;
    2816      649362 :   if (lgefint(p) == 3)
    2817             :   {
    2818      640731 :     ulong pp = uel(p, 2);
    2819      640731 :     r = Flx_to_ZX_inplace(Flx_rem(RgX_to_Flx(x, pp),
    2820             :                                   RgX_to_Flx(y, pp), pp));
    2821             :   }
    2822             :   else
    2823        8631 :     r = FpX_rem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2824      649362 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2825      451738 :   return gerepileupto(av, FpX_to_mod(r, p));
    2826             : }
    2827             : 
    2828             : static GEN
    2829       37631 : RgX_rem_QXQX(GEN x, GEN y, GEN T)
    2830             : {
    2831       37631 :   pari_sp av = avma;
    2832             :   GEN r;
    2833       37631 :   r = RgXQX_rem(RgX_liftred(x, T), RgX_liftred(y, T), T);
    2834       37631 :   return gerepilecopy(av, QXQX_to_mod_shallow(r, T));
    2835             : }
    2836             : static GEN
    2837          42 : RgX_rem_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    2838             : {
    2839          42 :   pari_sp av = avma;
    2840             :   GEN r;
    2841          42 :   GEN T = RgX_to_FpX(pol, p);
    2842          42 :   if (signe(T) == 0) pari_err_OP("%", x, y);
    2843          35 :   if (lgefint(p) == 3)
    2844             :   {
    2845          35 :     ulong pp = uel(p, 2);
    2846          35 :     GEN Tp = ZX_to_Flx(T, pp);
    2847          35 :     r = FlxX_to_ZXX(FlxqX_rem(RgX_to_FlxqX(x, Tp, pp),
    2848             :                               RgX_to_FlxqX(y, Tp, pp), Tp, pp));
    2849             :   }
    2850             :   else
    2851           0 :     r = FpXQX_rem(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), T, p);
    2852          35 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    2853          28 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    2854             : }
    2855             : 
    2856             : #define code(t1,t2) ((t1 << 6) | t2)
    2857             : static GEN
    2858    44853185 : RgX_mul_fast(GEN x, GEN y)
    2859             : {
    2860             :   GEN p, pol;
    2861             :   long pa;
    2862    44853185 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2863    44853185 :   switch(t)
    2864             :   {
    2865    33415078 :     case t_INT:    return ZX_mul(x,y);
    2866      985403 :     case t_FRAC:   return QX_mul(x,y);
    2867      101760 :     case t_FFELT:  return FFX_mul(x, y, pol);
    2868      494605 :     case t_INTMOD: return RgX_mul_FpX(x, y, p);
    2869             :     case code(t_POLMOD, t_INT):
    2870             :     case code(t_POLMOD, t_FRAC):
    2871      128077 :                    return RgX_mul_QXQX(x, y, pol);
    2872             :     case code(t_POLMOD, t_INTMOD):
    2873          77 :                    return RgX_mul_FpXQX(x, y, pol, p);
    2874     9728185 :     default:       return NULL;
    2875             :   }
    2876             : }
    2877             : static GEN
    2878      132196 : RgX_sqr_fast(GEN x)
    2879             : {
    2880             :   GEN p, pol;
    2881             :   long pa;
    2882      132196 :   long t = RgX_type(x,&p,&pol,&pa);
    2883      132196 :   switch(t)
    2884             :   {
    2885       86281 :     case t_INT:    return ZX_sqr(x);
    2886       32669 :     case t_FRAC:   return QX_sqr(x);
    2887        2261 :     case t_FFELT:  return FFX_sqr(x, pol);
    2888        1260 :     case t_INTMOD: return RgX_sqr_FpX(x, p);
    2889             :     case code(t_POLMOD, t_INT):
    2890             :     case code(t_POLMOD, t_FRAC):
    2891        3768 :                    return RgX_sqr_QXQX(x, pol);
    2892             :     case code(t_POLMOD, t_INTMOD):
    2893         154 :                    return RgX_sqr_FpXQX(x, pol, p);
    2894        5803 :     default:       return NULL;
    2895             :   }
    2896             : }
    2897             : 
    2898             : static GEN
    2899     5871324 : RgX_rem_fast(GEN x, GEN y)
    2900             : {
    2901             :   GEN p, pol;
    2902             :   long pa;
    2903     5871324 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    2904     5871324 :   switch(t)
    2905             :   {
    2906     3045866 :     case t_INT:    return ZX_is_monic(y) ? ZX_rem(x,y): NULL;
    2907     1223221 :     case t_FRAC:   return RgX_is_ZX(y) && ZX_is_monic(y) ? QX_ZX_rem(x,y): NULL;
    2908        1883 :     case t_FFELT:  return FFX_rem(x, y, pol);
    2909      649362 :     case t_INTMOD: return RgX_rem_FpX(x, y, p);
    2910             :     case code(t_POLMOD, t_INT):
    2911             :     case code(t_POLMOD, t_FRAC):
    2912       37631 :                    return RgX_rem_QXQX(x, y, pol);
    2913             :     case code(t_POLMOD, t_INTMOD):
    2914          42 :                    return RgX_rem_FpXQX(x, y, pol, p);
    2915      913319 :     default:       return NULL;
    2916             :   }
    2917             : }
    2918             : 
    2919             : #undef code
    2920             : 
    2921             : GEN
    2922    42337565 : RgX_mul(GEN x, GEN y)
    2923             : {
    2924    42337565 :   GEN z = RgX_mul_fast(x,y);
    2925    42337553 :   if (!z) z = RgX_mul_i(x,y);
    2926    42337553 :   return z;
    2927             : }
    2928             : 
    2929             : GEN
    2930      115557 : RgX_sqr(GEN x)
    2931             : {
    2932      115557 :   GEN z = RgX_sqr_fast(x);
    2933      115550 :   if (!z) z = RgX_sqr_i(x);
    2934      115550 :   return z;
    2935             : }
    2936             : 
    2937             : GEN
    2938     5871324 : RgX_rem(GEN x, GEN y)
    2939             : {
    2940     5871324 :   GEN z = RgX_rem_fast(x, y);
    2941     5871317 :   if (!z) z = RgX_divrem_i(x, y, ONLY_REM);
    2942     5871317 :   return z;
    2943             : }
    2944             : 
    2945             : GEN
    2946     2472032 : RgXn_mul(GEN f, GEN g, long n)
    2947             : {
    2948     2472032 :   pari_sp av = avma;
    2949     2472032 :   GEN h = RgX_mul_fast(f,g);
    2950     2472032 :   if (!h) return RgXn_mul2(f,g,n);
    2951      165236 :   if (degpol(h) < n) return h;
    2952       41440 :   return gerepilecopy(av, RgXn_red_shallow(h, n));
    2953             : }
    2954             : 
    2955             : GEN
    2956       16639 : RgXn_sqr(GEN f, long n)
    2957             : {
    2958       16639 :   pari_sp av = avma;
    2959       16639 :   GEN g = RgX_sqr_fast(f);
    2960       16639 :   if (!g) return RgXn_sqr2(f,n);
    2961       16401 :   if (degpol(g) < n) return g;
    2962        1071 :   return gerepilecopy(av, RgXn_red_shallow(g, n));
    2963             : }
    2964             : 
    2965             : /* (f*g) \/ x^n */
    2966             : GEN
    2967       43588 : RgX_mulhigh_i(GEN f, GEN g, long n)
    2968             : {
    2969       43588 :   GEN h = RgX_mul_fast(f,g);
    2970       43588 :   return h? RgX_shift_shallow(h, -n): RgX_mulhigh_i2(f,g,n);
    2971             : }
    2972             : 
    2973             : /* (f*g) \/ x^n */
    2974             : GEN
    2975           0 : RgX_sqrhigh_i(GEN f, long n)
    2976             : {
    2977           0 :   GEN h = RgX_sqr_fast(f);
    2978           0 :   return h? RgX_shift_shallow(h, -n): RgX_sqrhigh_i2(f,n);
    2979             : }

Generated by: LCOV version 1.13