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.18.1 lcov report (development 30005-fc14bb602a) Lines: 1628 1784 91.3 %
Date: 2025-02-18 09:22:46 Functions: 204 222 91.9 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation; either version 2 of the License, or (at your option) any later
       8             : version. It is distributed in the hope that it will be useful, but WITHOUT
       9             : ANY WARRANTY WHATSOEVER.
      10             : 
      11             : Check the License for details. You should have received a copy of it, along
      12             : with the package; see the file 'COPYING'. If not, write to the Free Software
      13             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      14             : 
      15             : #include "pari.h"
      16             : #include "paripriv.h"
      17             : 
      18             : #define DEBUGLEVEL DEBUGLEVEL_pol
      19             : 
      20             : /*******************************************************************/
      21             : /*                                                                 */
      22             : /*                         GENERIC                                 */
      23             : /*                                                                 */
      24             : /*******************************************************************/
      25             : 
      26             : /* Return optimal parameter l for the evaluation of n/m polynomials of degree d
      27             :    Fractional values can be used if the evaluations are done with different
      28             :    accuracies, and thus have different weights.
      29             :  */
      30             : long
      31    19460170 : brent_kung_optpow(long d, long n, long m)
      32             : {
      33             :   long p, r;
      34    19460170 :   long pold=1, rold=n*(d-1);
      35    91235846 :   for(p=2; p<=d; p++)
      36             :   {
      37    71775676 :     r = m*(p-1) + n*((d-1)/p);
      38    71775676 :     if (r<rold) { pold=p; rold=r; }
      39             :   }
      40    19460170 :   return pold;
      41             : }
      42             : 
      43             : static GEN
      44    14563704 : gen_RgXQ_eval_powers(GEN P, GEN V, long a, long n, void *E, const struct bb_algebra *ff,
      45             :                                            GEN cmul(void *E, GEN P, long a, GEN x))
      46             : {
      47    14563704 :   pari_sp av = avma;
      48             :   long i;
      49    14563704 :   GEN z = cmul(E,P,a,ff->one(E));
      50    14537959 :   if (!z) z = gen_0;
      51    76359433 :   for (i=1; i<=n; i++)
      52             :   {
      53    61826782 :     GEN t = cmul(E,P,a+i,gel(V,i+1));
      54    61853249 :     if (t) {
      55    46932490 :       z = ff->add(E, z, t);
      56    46896871 :       if (gc_needed(av,2)) z = gerepileupto(av, z);
      57             :     }
      58             :   }
      59    14532651 :   return ff->red(E,z);
      60             : }
      61             : 
      62             : /* Brent & Kung
      63             :  * (Fast algorithms for manipulating formal power series, JACM 25:581-595, 1978)
      64             :  *
      65             :  * V as output by FpXQ_powers(x,l,T,p). For optimal performance, l is as given
      66             :  * by brent_kung_optpow */
      67             : GEN
      68    12470107 : gen_bkeval_powers(GEN P, long d, GEN V, void *E, const struct bb_algebra *ff,
      69             :                                      GEN cmul(void *E, GEN P, long a, GEN x))
      70             : {
      71    12470107 :   pari_sp av = avma;
      72    12470107 :   long l = lg(V)-1;
      73             :   GEN z, u;
      74             : 
      75    12470107 :   if (d < 0) return ff->zero(E);
      76    11862481 :   if (d < l) return gerepileupto(av, gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul));
      77     1158040 :   if (l<2) pari_err_DOMAIN("gen_RgX_bkeval_powers", "#powers", "<",gen_2,V);
      78     1158040 :   if (DEBUGLEVEL>=8)
      79             :   {
      80           0 :     long cnt = 1 + (d - l) / (l-1);
      81           0 :     err_printf("RgX_RgXQV_eval(%ld/%ld): %ld RgXQ_mul\n", d, l-1, cnt);
      82             :   }
      83     1158040 :   d -= l;
      84     1158040 :   z = gen_RgXQ_eval_powers(P,V,d+1,l-1,E,ff,cmul);
      85     2700149 :   while (d >= l-1)
      86             :   {
      87     1540964 :     d -= l-1;
      88     1540964 :     u = gen_RgXQ_eval_powers(P,V,d+1,l-2,E,ff,cmul);
      89     1540926 :     z = ff->add(E,u, ff->mul(E,z,gel(V,l)));
      90     1540959 :     if (gc_needed(av,2))
      91          91 :       z = gerepileupto(av, z);
      92             :   }
      93     1159185 :   u = gen_RgXQ_eval_powers(P,V,0,d,E,ff,cmul);
      94     1159196 :   z = ff->add(E,u, ff->mul(E,z,gel(V,d+2)));
      95     1159189 :   return gerepileupto(av, ff->red(E,z));
      96             : }
      97             : 
      98             : GEN
      99      874097 : gen_bkeval(GEN Q, long d, GEN x, int use_sqr, void *E, const struct bb_algebra *ff,
     100             :                                       GEN cmul(void *E, GEN P, long a, GEN x))
     101             : {
     102      874097 :   pari_sp av = avma;
     103             :   GEN z, V;
     104             :   long rtd;
     105      874097 :   if (d < 0) return ff->zero(E);
     106      872802 :   rtd = (long) sqrt((double)d);
     107      872802 :   V = gen_powers(x,rtd,use_sqr,E,ff->sqr,ff->mul,ff->one);
     108      872811 :   z = gen_bkeval_powers(Q, d, V, E, ff, cmul);
     109      872806 :   return gerepileupto(av, z);
     110             : }
     111             : 
     112             : static GEN
     113     2027284 : _gen_nored(void *E, GEN x) { (void)E; return x; }
     114             : static GEN
     115    19235852 : _gen_add(void *E, GEN x, GEN y) { (void)E; return gadd(x, y); }
     116             : static GEN
     117           0 : _gen_sub(void *E, GEN x, GEN y) { (void)E; return gsub(x, y); }
     118             : static GEN
     119     1839431 : _gen_mul(void *E, GEN x, GEN y) { (void)E; return gmul(x, y); }
     120             : static GEN
     121      591635 : _gen_sqr(void *E, GEN x) { (void)E; return gsqr(x); }
     122             : static GEN
     123     2068456 : _gen_one(void *E) { (void)E; return gen_1; }
     124             : static GEN
     125       24184 : _gen_zero(void *E) { (void)E; return gen_0; }
     126             : 
     127             : static struct bb_algebra Rg_algebra = { _gen_nored, _gen_add, _gen_sub,
     128             :               _gen_mul, _gen_sqr,_gen_one,_gen_zero };
     129             : 
     130             : static GEN
     131      512818 : _gen_cmul(void *E, GEN P, long a, GEN x)
     132      512818 : {(void)E; return gmul(gel(P,a+2), x);}
     133             : 
     134             : GEN
     135      166768 : RgX_RgV_eval(GEN Q, GEN x)
     136             : {
     137      166768 :   return gen_bkeval_powers(Q, degpol(Q), x, NULL, &Rg_algebra, _gen_cmul);
     138             : }
     139             : 
     140             : GEN
     141           0 : RgX_Rg_eval_bk(GEN Q, GEN x)
     142             : {
     143           0 :   return gen_bkeval(Q, degpol(Q), x, 1, NULL, &Rg_algebra, _gen_cmul);
     144             : }
     145             : 
     146             : GEN
     147        2947 : RgXV_RgV_eval(GEN Q, GEN x)
     148             : {
     149        2947 :   long i, l = lg(Q), vQ = gvar(Q);
     150        2947 :   GEN v = cgetg(l, t_VEC);
     151      248311 :   for (i = 1; i < l; i++)
     152             :   {
     153      245364 :     GEN Qi = gel(Q, i);
     154      245364 :     gel(v, i) = typ(Qi)==t_POL && varn(Qi)==vQ? RgX_RgV_eval(Qi, x): gcopy(Qi);
     155             :   }
     156        2947 :   return v;
     157             : }
     158             : 
     159             : GEN
     160      671155 : RgX_homogenous_evalpow(GEN P, GEN A, GEN B)
     161             : {
     162      671155 :   pari_sp av = avma, btop;
     163      671155 :   long i, d = degpol(P), o;
     164             :   GEN s;
     165      671154 :   if (signe(P)==0) return pol_0(varn(P));
     166      671154 :   s = gel(P, d+2);
     167      671154 :   if (d == 0) return gcopy(s);
     168      667878 :   o = RgX_deflate_order(P);
     169      667900 :   if (o > 1) A = gpowgs(A, o);
     170      667906 :   btop = avma;
     171     2318560 :   for (i = d-o; i >= 0; i-=o)
     172             :   {
     173     1650660 :     s = gadd(gmul(s, A), gmul(gel(B,d+1-i), gel(P,i+2)));
     174     1650654 :     if (gc_needed(btop,1))
     175             :     {
     176          13 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_homogenous_eval(%ld)",i);
     177          13 :       s = gerepileupto(btop, s);
     178             :     }
     179             :   }
     180      667900 :   return gerepileupto(av, s);
     181             : }
     182             : 
     183             : GEN
     184        1652 : QXQX_homogenous_evalpow(GEN P, GEN A, GEN B, GEN T)
     185             : {
     186        1652 :   pari_sp av = avma;
     187        1652 :   long i, d = degpol(P), v = varn(A);
     188             :   GEN s;
     189        1652 :   if (signe(P)==0) return pol_0(v);
     190        1652 :   if (d == 0) return scalarpol(gel(P, d+2), v);
     191        1232 :   s = scalarpol_shallow(gel(P, d+2), v);
     192        4963 :   for (i = d-1; i >= 0; i--)
     193             :   {
     194        3731 :     GEN c = gel(P,i+2), b = gel(B,d+1-i);
     195        3731 :     s = RgX_add(QXQX_mul(s, A, T), typ(c)==t_POL ? QXQX_QXQ_mul(b, c, T): gmul(b, c));
     196        3731 :     if (gc_needed(av,1))
     197             :     {
     198           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"QXQX_homogenous_eval(%ld)",i);
     199           0 :       s = gerepileupto(av, s);
     200             :     }
     201             :   }
     202        1232 :   return gerepileupto(av, s);
     203             : }
     204             : 
     205             : const struct bb_algebra *
     206      281039 : get_Rg_algebra(void)
     207             : {
     208      281039 :   return &Rg_algebra;
     209             : }
     210             : 
     211             : static struct bb_ring Rg_ring = {  _gen_add, _gen_mul, _gen_sqr };
     212             : 
     213             : static GEN
     214       11333 : _RgX_divrem(void *E, GEN x, GEN y, GEN *r)
     215             : {
     216             :   (void) E;
     217       11333 :   return RgX_divrem(x, y, r);
     218             : }
     219             : 
     220             : GEN
     221        3157 : RgX_digits(GEN x, GEN T)
     222             : {
     223        3157 :   long d = degpol(T), n = (lgpol(x)+d-1)/d;
     224        3157 :   return gen_digits(x,T,n,NULL, &Rg_ring, _RgX_divrem);
     225             : }
     226             : 
     227             : /*******************************************************************/
     228             : /*                                                                 */
     229             : /*                         RgX                                     */
     230             : /*                                                                 */
     231             : /*******************************************************************/
     232             : 
     233             : long
     234    24847703 : RgX_equal(GEN x, GEN y)
     235             : {
     236    24847703 :   long i = lg(x);
     237             : 
     238    24847703 :   if (i != lg(y)) return 0;
     239   106959197 :   for (i--; i > 1; i--)
     240    82439745 :     if (!gequal(gel(x,i),gel(y,i))) return 0;
     241    24519452 :   return 1;
     242             : }
     243             : 
     244             : /* Returns 1 in the base ring over which x is defined */
     245             : /* HACK: this also works for t_SER */
     246             : GEN
     247   156430459 : Rg_get_1(GEN x)
     248             : {
     249             :   GEN p, T;
     250   156430459 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     251   156430459 :   if (RgX_type_is_composite(tx))
     252    11648957 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     253   156430459 :   switch(tx)
     254             :   {
     255      427217 :     case t_INTMOD: retmkintmod(is_pm1(p)? gen_0: gen_1, icopy(p));
     256         945 :     case t_PADIC: return cvtop(gen_1, p, lx);
     257        4998 :     case t_FFELT: return FF_1(T);
     258   155997299 :     default: return gen_1;
     259             :   }
     260             : }
     261             : /* Returns 0 in the base ring over which x is defined */
     262             : /* HACK: this also works for t_SER */
     263             : GEN
     264     6994197 : Rg_get_0(GEN x)
     265             : {
     266             :   GEN p, T;
     267     6994197 :   long i, lx, tx = Rg_type(x, &p, &T, &lx);
     268     6994197 :   if (RgX_type_is_composite(tx))
     269       60732 :     RgX_type_decode(tx, &i /*junk*/, &tx);
     270     6994197 :   switch(tx)
     271             :   {
     272         497 :     case t_INTMOD: retmkintmod(gen_0, icopy(p));
     273          42 :     case t_PADIC: return zeropadic(p, lx);
     274         231 :     case t_FFELT: return FF_zero(T);
     275     6993427 :     default: return gen_0;
     276             :   }
     277             : }
     278             : 
     279             : GEN
     280        6902 : QX_ZXQV_eval(GEN P, GEN V, GEN dV)
     281             : {
     282        6902 :   long i, n = degpol(P);
     283             :   GEN z, dz, dP;
     284        6902 :   if (n < 0) return gen_0;
     285        6902 :   P = Q_remove_denom(P, &dP);
     286        6902 :   z = gel(P,2); if (n == 0) return icopy(z);
     287        3906 :   if (dV) z = mulii(dV, z); /* V[1] = dV */
     288        3906 :   z = ZX_Z_add_shallow(ZX_Z_mul(gel(V,2),gel(P,3)), z);
     289        7525 :   for (i=2; i<=n; i++) z = ZX_add(ZX_Z_mul(gel(V,i+1),gel(P,2+i)), z);
     290        3906 :   dz = mul_denom(dP, dV);
     291        3906 :   return dz? RgX_Rg_div(z, dz): z;
     292             : }
     293             : 
     294             : /* Return P(h * x), not memory clean */
     295             : GEN
     296       50920 : RgX_unscale(GEN P, GEN h)
     297             : {
     298       50920 :   long i, l = lg(P);
     299       50920 :   GEN hi = gen_1, Q = cgetg(l, t_POL);
     300       50920 :   Q[1] = P[1];
     301       50920 :   if (l == 2) return Q;
     302       39930 :   gel(Q,2) = gcopy(gel(P,2));
     303       94070 :   for (i=3; i<l; i++)
     304             :   {
     305       54140 :     hi = gmul(hi,h);
     306       54140 :     gel(Q,i) = gmul(gel(P,i), hi);
     307             :   }
     308       39930 :   return Q;
     309             : }
     310             : /* P a ZX, Return P(h * x), not memory clean; optimize for h = -1 */
     311             : GEN
     312     1473362 : ZX_z_unscale(GEN P, long h)
     313             : {
     314     1473362 :   long i, l = lg(P);
     315     1473362 :   GEN Q = cgetg(l, t_POL);
     316     1473362 :   Q[1] = P[1];
     317     1473362 :   if (l == 2) return Q;
     318     1384816 :   gel(Q,2) = gel(P,2);
     319     1384816 :   if (l == 3) return Q;
     320     1358510 :   if (h == -1)
     321      237548 :     for (i = 3; i < l; i++)
     322             :     {
     323      196827 :       gel(Q,i) = negi(gel(P,i));
     324      196827 :       if (++i == l) break;
     325      145799 :       gel(Q,i) = gel(P,i);
     326             :     }
     327             :   else
     328             :   {
     329             :     GEN hi;
     330     1266761 :     gel(Q,3) = mulis(gel(P,3), h);
     331     1266761 :     hi = sqrs(h);
     332     5450372 :     for (i = 4; i < l; i++)
     333             :     {
     334     4183611 :       gel(Q,i) = mulii(gel(P,i), hi);
     335     4183611 :       if (i != l-1) hi = mulis(hi,h);
     336             :     }
     337             :   }
     338     1358510 :   return Q;
     339             : }
     340             : /* P a ZX, h a t_INT. Return P(h * x), not memory clean; optimize for h = -1 */
     341             : GEN
     342     1019154 : ZX_unscale(GEN P, GEN h)
     343             : {
     344             :   long i, l;
     345             :   GEN Q, hi;
     346     1019154 :   i = itos_or_0(h); if (i) return ZX_z_unscale(P, i);
     347         881 :   l = lg(P); Q = cgetg(l, t_POL);
     348         881 :   Q[1] = P[1];
     349         881 :   if (l == 2) return Q;
     350         881 :   gel(Q,2) = gel(P,2);
     351         881 :   if (l == 3) return Q;
     352         881 :   hi = h;
     353         881 :   gel(Q,3) = mulii(gel(P,3), hi);
     354        2741 :   for (i = 4; i < l; i++)
     355             :   {
     356        1860 :     hi = mulii(hi,h);
     357        1860 :     gel(Q,i) = mulii(gel(P,i), hi);
     358             :   }
     359         881 :   return Q;
     360             : }
     361             : /* P a ZX. Return P(x << n), not memory clean */
     362             : GEN
     363      927756 : ZX_unscale2n(GEN P, long n)
     364             : {
     365      927756 :   long i, ni = n, l = lg(P);
     366      927756 :   GEN Q = cgetg(l, t_POL);
     367      927756 :   Q[1] = P[1];
     368      927756 :   if (l == 2) return Q;
     369      927756 :   gel(Q,2) = gel(P,2);
     370      927756 :   if (l == 3) return Q;
     371      927756 :   gel(Q,3) = shifti(gel(P,3), ni);
     372     3172065 :   for (i=4; i<l; i++)
     373             :   {
     374     2244368 :     ni += n;
     375     2244368 :     gel(Q,i) = shifti(gel(P,i), ni);
     376             :   }
     377      927697 :   return Q;
     378             : }
     379             : /* P(h*X) / h, assuming h | P(0), i.e. the result is a ZX */
     380             : GEN
     381       12997 : ZX_unscale_div(GEN P, GEN h)
     382             : {
     383       12997 :   long i, l = lg(P);
     384       12997 :   GEN hi, Q = cgetg(l, t_POL);
     385       12997 :   Q[1] = P[1];
     386       12997 :   if (l == 2) return Q;
     387       12997 :   gel(Q,2) = diviiexact(gel(P,2), h);
     388       12997 :   if (l == 3) return Q;
     389       12997 :   gel(Q,3) = gel(P,3);
     390       12997 :   if (l == 4) return Q;
     391       12997 :   hi = h;
     392       12997 :   gel(Q,4) = mulii(gel(P,4), hi);
     393       64121 :   for (i=5; i<l; i++)
     394             :   {
     395       51124 :     hi = mulii(hi,h);
     396       51124 :     gel(Q,i) = mulii(gel(P,i), hi);
     397             :   }
     398       12997 :   return Q;
     399             : }
     400             : /* P(h*X) / h^k, assuming the result is a ZX */
     401             : GEN
     402        1393 : ZX_unscale_divpow(GEN P, GEN h, long k)
     403             : {
     404        1393 :   long i, j, l = lg(P);
     405        1393 :   GEN H, Q = cgetg(l, t_POL);
     406        1393 :   Q[1] = P[1]; if (l == 2) return Q;
     407        1393 :   H = gpowers(h, maxss(k, l - 3 - k));
     408        5572 :   for (i = 2, j = k+1; j > 1 && i < l; i++)
     409        4179 :     gel(Q, i) = diviiexact(gel(P, i), gel(H, j--));
     410        1393 :   if (i == l) return Q;
     411        1393 :   gel(Q, i) = gel(P, i); i++;
     412        5082 :   for (j = 2; i < l; i++) gel(Q, i) = mulii(gel(P, i), gel(H, j++));
     413        1393 :   return Q;
     414             : }
     415             : 
     416             : GEN
     417        6251 : RgXV_unscale(GEN x, GEN h)
     418             : {
     419        6251 :   if (isint1(h)) return gcopy(x);
     420       18320 :   pari_APPLY_same(RgX_unscale(gel(x,i), h));
     421             : }
     422             : 
     423             : /* Return h^degpol(P) P(x / h), not memory clean */
     424             : GEN
     425     4343829 : RgX_rescale(GEN P, GEN h)
     426             : {
     427     4343829 :   long i, l = lg(P);
     428     4343829 :   GEN Q = cgetg(l,t_POL), hi = h;
     429     4343833 :   gel(Q,l-1) = gel(P,l-1);
     430    11431920 :   for (i=l-2; i>=2; i--)
     431             :   {
     432    11429288 :     gel(Q,i) = gmul(gel(P,i), hi);
     433    11429141 :     if (i == 2) break;
     434     7087756 :     hi = gmul(hi,h);
     435             :   }
     436     4344017 :   Q[1] = P[1]; return Q;
     437             : }
     438             : 
     439             : GEN
     440        2401 : RgXV_rescale(GEN x, GEN h)
     441             : {
     442        2401 :   if (isint1(h)) return RgX_copy(x);
     443       16086 :   pari_APPLY_same(RgX_rescale(gel(x,i), h));
     444             : }
     445             : 
     446             : /* A(X^d) --> A(X) */
     447             : GEN
     448      987833 : RgX_deflate(GEN x0, long d)
     449             : {
     450             :   GEN z, y, x;
     451      987833 :   long i,id, dy, dx = degpol(x0);
     452      987833 :   if (d == 1 || dx <= 0) return leafcopy(x0);
     453      376936 :   dy = dx/d;
     454      376936 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     455      376936 :   z = y + 2;
     456      376936 :   x = x0+ 2;
     457     1516593 :   for (i=id=0; i<=dy; i++,id+=d) gel(z,i) = gel(x,id);
     458      376936 :   return y;
     459             : }
     460             : 
     461             : GEN
     462        1260 : RgX_homogenize_deg(GEN P, long d, long v)
     463             : {
     464             :   long i, l;
     465        1260 :   GEN Q = cgetg_copy(P, &l);
     466        1260 :   Q[1] = P[1];
     467        4018 :   for (i = 2; i < l; i++) gel(Q,i) = monomial(gel(P,i), d--, v);
     468        1260 :   return Q;
     469             : }
     470             : 
     471             : GEN
     472       18732 : RgX_homogenize(GEN P, long v)
     473             : {
     474             :   long i, l, d;
     475       18732 :   GEN Q = cgetg_copy(P, &l);
     476       18732 :   Q[1] = P[1]; d = l-3;
     477      165998 :   for (i = 2; i < l; i++) gel(Q,i) = monomial(gel(P,i), d--, v);
     478       18732 :   return Q;
     479             : }
     480             : 
     481             : /* F a t_RFRAC */
     482             : long
     483         105 : rfrac_deflate_order(GEN F)
     484             : {
     485         105 :   GEN N = gel(F,1), D = gel(F,2);
     486         105 :   long m = (degpol(D) <= 0)? 0: RgX_deflate_order(D);
     487         105 :   if (m == 1) return 1;
     488          42 :   if (typ(N) == t_POL && varn(N) == varn(D))
     489          28 :     m = cgcd(m, RgX_deflate_order(N));
     490          42 :   return m;
     491             : }
     492             : /* F a t_RFRAC */
     493             : GEN
     494         105 : rfrac_deflate_max(GEN F, long *m)
     495             : {
     496         105 :   *m = rfrac_deflate_order(F);
     497         105 :   return rfrac_deflate(F, *m);
     498             : }
     499             : /* F a t_RFRAC */
     500             : GEN
     501         105 : rfrac_deflate(GEN F, long m)
     502             : {
     503         105 :   GEN N = gel(F,1), D = gel(F,2);
     504         105 :   if (m == 1) return F;
     505          42 :   if (typ(N) == t_POL && varn(N) == varn(D)) N = RgX_deflate(N, m);
     506          42 :   D = RgX_deflate(D, m); return mkrfrac(N, D);
     507             : }
     508             : 
     509             : /* return x0(X^d) */
     510             : GEN
     511      798877 : RgX_inflate(GEN x0, long d)
     512             : {
     513      798877 :   long i, id, dy, dx = degpol(x0);
     514      798877 :   GEN x = x0 + 2, z, y;
     515      798877 :   if (dx <= 0) return leafcopy(x0);
     516      791261 :   dy = dx*d;
     517      791261 :   y = cgetg(dy+3, t_POL); y[1] = x0[1];
     518      791261 :   z = y + 2;
     519    27857323 :   for (i=0; i<=dy; i++) gel(z,i) = gen_0;
     520    11411892 :   for (i=id=0; i<=dx; i++,id+=d) gel(z,id) = gel(x,i);
     521      791261 :   return y;
     522             : }
     523             : 
     524             : /* return P(X + c) using destructive Horner, optimize for c = 1,-1 */
     525             : static GEN
     526     5699886 : RgX_translate_basecase(GEN P, GEN c)
     527             : {
     528     5699886 :   pari_sp av = avma;
     529             :   GEN Q, R;
     530             :   long i, k, n;
     531             : 
     532     5699886 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     533     5699864 :   Q = leafcopy(P);
     534     5699861 :   R = Q+2; n = degpol(P);
     535     5699858 :   if (isint1(c))
     536             :   {
     537        7595 :     for (i=1; i<=n; i++)
     538             :     {
     539       20265 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gel(R,k+1));
     540        5390 :       if (gc_needed(av,2))
     541             :       {
     542           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(1), i = %ld/%ld", i,n);
     543           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     544             :       }
     545             :     }
     546             :   }
     547     5697652 :   else if (isintm1(c))
     548             :   {
     549       15974 :     for (i=1; i<=n; i++)
     550             :     {
     551       49630 :       for (k=n-i; k<n; k++) gel(R,k) = gsub(gel(R,k), gel(R,k+1));
     552       12299 :       if (gc_needed(av,2))
     553             :       {
     554           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate(-1), i = %ld/%ld", i,n);
     555           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     556             :       }
     557             :     }
     558             :   }
     559             :   else
     560             :   {
     561    19427406 :     for (i=1; i<=n; i++)
     562             :     {
     563    45467434 :       for (k=n-i; k<n; k++) gel(R,k) = gadd(gel(R,k), gmul(c, gel(R,k+1)));
     564    13733429 :       if (gc_needed(av,2))
     565             :       {
     566           0 :         if(DEBUGMEM>1) pari_warn(warnmem,"RgX_translate, i = %ld/%ld", i,n);
     567           0 :         Q = gerepilecopy(av, Q); R = Q+2;
     568             :       }
     569             :     }
     570             :   }
     571     5699567 :   return gerepilecopy(av, Q);
     572             : }
     573             : 
     574             : static GEN
     575      142269 : zero_FpX_mod(GEN p, long v)
     576             : {
     577      142269 :   GEN r = cgetg(3,t_POL);
     578      142269 :   r[1] = evalvarn(v);
     579      142269 :   gel(r,2) = mkintmod(gen_0, icopy(p));
     580      142269 :   return r;
     581             : }
     582             : 
     583             : static GEN
     584           0 : RgX_translate_FpX(GEN P, GEN c, GEN p)
     585             : {
     586           0 :   pari_sp av = avma;
     587             :   GEN r;
     588             : #if 0
     589             :   if (lgefint(p) == 3)
     590             :   {
     591             :     ulong pp = uel(p, 2);
     592             :     r = Flx_to_ZX_inplace(Flx_translate(RgX_to_Flx(x, pp), Rg_to_Fl(c, pp), pp));
     593             :   }
     594             :   else
     595             : #endif
     596           0 :     r = FpX_translate(RgX_to_FpX(P, p), Rg_to_Fp(c, p), p);
     597           0 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(P)); }
     598           0 :   return gerepileupto(av, FpX_to_mod(r, p));
     599             : }
     600             : 
     601             : static GEN
     602     5982199 : RgX_translate_fast(GEN P, GEN c)
     603             : {
     604             :   GEN p, pol;
     605             :   long pa;
     606     5982199 :   long t = RgX_Rg_type(P, c, &p,&pol,&pa);
     607     5982198 :   switch(t)
     608             :   {
     609      282502 :     case t_INT:    return ZX_translate(P, c);
     610           0 :     case t_INTMOD: return RgX_translate_FpX(P, c, p);
     611     5699698 :     default:       return NULL;
     612             :   }
     613             : }
     614             : 
     615             : static GEN
     616     5700075 : RgX_translate_i(GEN P, GEN c)
     617             : {
     618     5700075 :   pari_sp av = avma;
     619             :   long n;
     620     5700075 :   n = degpol(P);
     621     5700074 :   if (n < 40)
     622     5699885 :     return RgX_translate_basecase(P, c);
     623             :   else
     624             :   {
     625         189 :     long d = n >> 1;
     626         189 :     GEN Q = RgX_translate_i(RgX_shift_shallow(P, -d), c);
     627         189 :     GEN R = RgX_translate_i(RgXn_red_shallow(P, d), c);
     628         189 :     GEN S = gpowgs(deg1pol_shallow(gen_1, c, varn(P)), d);
     629         189 :     return gerepileupto(av, RgX_add(RgX_mul(Q, S), R));
     630             :   }
     631             : }
     632             : 
     633             : GEN
     634     5982201 : RgX_translate(GEN P, GEN c)
     635             : {
     636     5982201 :   GEN R = RgX_translate_fast(P, c);
     637     5982202 :   return R ? R: RgX_translate_i(P,c);
     638             : }
     639             : /* P(ax + b) */
     640             : GEN
     641       29666 : RgX_affine(GEN P, GEN a, GEN b)
     642             : {
     643       29666 :   if (!gequal0(b)) P = RgX_translate(P, b);
     644       29666 :   return RgX_unscale(P, a);
     645             : }
     646             : 
     647             : /* return lift( P(X + c) ) using Horner, c in R[y]/(T) */
     648             : GEN
     649       32464 : RgXQX_translate(GEN P, GEN c, GEN T)
     650             : {
     651       32464 :   pari_sp av = avma;
     652             :   GEN Q, R;
     653             :   long i, k, n;
     654             : 
     655       32464 :   if (!signe(P) || gequal0(c)) return RgX_copy(P);
     656       32121 :   Q = leafcopy(P);
     657       32121 :   R = Q+2; n = degpol(P);
     658      103762 :   for (i=1; i<=n; i++)
     659             :   {
     660      301670 :     for (k=n-i; k<n; k++)
     661             :     {
     662      230029 :       pari_sp av2 = avma;
     663      230029 :       gel(R,k) = gerepileupto(av2,
     664      230029 :                    RgX_rem(gadd(gel(R,k), gmul(c, gel(R,k+1))), T));
     665             :     }
     666       71641 :     if (gc_needed(av,2))
     667             :     {
     668           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXQX_translate, i = %ld/%ld", i,n);
     669           0 :       Q = gerepilecopy(av, Q); R = Q+2;
     670             :     }
     671             :   }
     672       32121 :   return gerepilecopy(av, Q);
     673             : }
     674             : 
     675             : /********************************************************************/
     676             : /**                                                                **/
     677             : /**                          CONVERSIONS                           **/
     678             : /**                       (not memory clean)                       **/
     679             : /**                                                                **/
     680             : /********************************************************************/
     681             : /* to INT / FRAC / (POLMOD mod T), not memory clean because T not copied,
     682             :  * but everything else is */
     683             : static GEN
     684      167790 : QXQ_to_mod(GEN x, GEN T)
     685             : {
     686             :   long d;
     687      167790 :   switch(typ(x))
     688             :   {
     689       62815 :     case t_INT:  return icopy(x);
     690        2079 :     case t_FRAC: return gcopy(x);
     691      102896 :     case t_POL:
     692      102896 :       d = degpol(x);
     693      102896 :       if (d < 0) return gen_0;
     694      100635 :       if (d == 0) return gcopy(gel(x,2));
     695       98482 :       return mkpolmod(RgX_copy(x), T);
     696           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     697             :              return NULL;/* LCOV_EXCL_LINE */
     698             :   }
     699             : }
     700             : /* pure shallow version */
     701             : GEN
     702      839020 : QXQ_to_mod_shallow(GEN x, GEN T)
     703             : {
     704             :   long d;
     705      839020 :   switch(typ(x))
     706             :   {
     707      542454 :     case t_INT:
     708      542454 :     case t_FRAC: return x;
     709      296566 :     case t_POL:
     710      296566 :       d = degpol(x);
     711      296566 :       if (d < 0) return gen_0;
     712      248992 :       if (d == 0) return gel(x,2);
     713      232033 :       return mkpolmod(x, T);
     714           0 :     default: pari_err_TYPE("QXQ_to_mod",x);
     715             :              return NULL;/* LCOV_EXCL_LINE */
     716             :   }
     717             : }
     718             : /* T a ZX, z lifted from (Q[Y]/(T(Y)))[X], apply QXQ_to_mod to all coeffs.
     719             :  * Not memory clean because T not copied, but everything else is */
     720             : static GEN
     721       37492 : QXQX_to_mod(GEN z, GEN T)
     722             : {
     723       37492 :   long i,l = lg(z);
     724       37492 :   GEN x = cgetg(l,t_POL);
     725      187033 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod(gel(z,i), T);
     726       37492 :   x[1] = z[1]; return normalizepol_lg(x,l);
     727             : }
     728             : /* pure shallow version */
     729             : GEN
     730      202707 : QXQX_to_mod_shallow(GEN z, GEN T)
     731             : {
     732      202707 :   long i,l = lg(z);
     733      202707 :   GEN x = cgetg(l,t_POL);
     734      958539 :   for (i=2; i<l; i++) gel(x,i) = QXQ_to_mod_shallow(gel(z,i), T);
     735      202707 :   x[1] = z[1]; return normalizepol_lg(x,l);
     736             : }
     737             : /* Apply QXQX_to_mod to all entries. Memory-clean ! */
     738             : GEN
     739       12656 : QXQXV_to_mod(GEN V, GEN T)
     740             : {
     741       12656 :   long i, l = lg(V);
     742       12656 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     743       50148 :   for (i=1;i<l; i++) gel(z,i) = QXQX_to_mod(gel(V,i), T);
     744       12656 :   return z;
     745             : }
     746             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     747             : GEN
     748       18370 : QXQV_to_mod(GEN V, GEN T)
     749             : {
     750       18370 :   long i, l = lg(V);
     751       18370 :   GEN z = cgetg(l, t_VEC); T = ZX_copy(T);
     752       36619 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod(gel(V,i), T);
     753       18370 :   return z;
     754             : }
     755             : 
     756             : /* Apply QXQ_to_mod to all entries. Memory-clean ! */
     757             : GEN
     758       14644 : QXQC_to_mod_shallow(GEN V, GEN T)
     759             : {
     760       14644 :   long i, l = lg(V);
     761       14644 :   GEN z = cgetg(l, t_COL);
     762       97832 :   for (i=1;i<l; i++) gel(z,i) = QXQ_to_mod_shallow(gel(V,i), T);
     763       14644 :   return z;
     764             : }
     765             : 
     766             : GEN
     767        6622 : QXQM_to_mod_shallow(GEN V, GEN T)
     768             : {
     769        6622 :   long i, l = lg(V);
     770        6622 :   GEN z = cgetg(l, t_MAT);
     771       21266 :   for (i=1; i<l; i++) gel(z,i) = QXQC_to_mod_shallow(gel(V,i), T);
     772        6622 :   return z;
     773             : }
     774             : 
     775             : GEN
     776     6167231 : RgX_renormalize_lg(GEN x, long lx)
     777             : {
     778             :   long i;
     779     9191801 :   for (i = lx-1; i>1; i--)
     780     8739575 :     if (! gequal0(gel(x,i))) break; /* _not_ isexactzero */
     781     6167231 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
     782     6167231 :   setlg(x, i+1); setsigne(x, i != 1); return x;
     783             : }
     784             : 
     785             : GEN
     786     1492872 : RgV_to_RgX(GEN x, long v)
     787             : {
     788     1492872 :   long i, k = lg(x);
     789             :   GEN p;
     790             : 
     791     4153856 :   while (--k && gequal0(gel(x,k)));
     792     1492875 :   if (!k) return pol_0(v);
     793     1485597 :   i = k+2; p = cgetg(i,t_POL);
     794     1485593 :   p[1] = evalsigne(1) | evalvarn(v);
     795    10231153 :   x--; for (k=2; k<i; k++) gel(p,k) = gel(x,k);
     796     1485593 :   return p;
     797             : }
     798             : GEN
     799      188283 : RgV_to_RgX_reverse(GEN x, long v)
     800             : {
     801      188283 :   long j, k, l = lg(x);
     802             :   GEN p;
     803             : 
     804      189851 :   for (k = 1; k < l; k++)
     805      189851 :     if (!gequal0(gel(x,k))) break;
     806      188283 :   if (k == l) return pol_0(v);
     807      188283 :   k -= 1;
     808      188283 :   l -= k;
     809      188283 :   x += k;
     810      188283 :   p = cgetg(l+1,t_POL);
     811      188283 :   p[1] = evalsigne(1) | evalvarn(v);
     812      985038 :   for (j=2, k=l; j<=l; j++) gel(p,j) = gel(x,--k);
     813      188283 :   return p;
     814             : }
     815             : 
     816             : /* return the (N-dimensional) vector of coeffs of p */
     817             : GEN
     818    13408527 : RgX_to_RgC(GEN x, long N)
     819             : {
     820             :   long i, l;
     821             :   GEN z;
     822    13408527 :   l = lg(x)-1; x++;
     823    13408527 :   if (l > N+1) l = N+1; /* truncate higher degree terms */
     824    13408527 :   z = cgetg(N+1,t_COL);
     825    80852118 :   for (i=1; i<l ; i++) gel(z,i) = gel(x,i);
     826    24467338 :   for (   ; i<=N; i++) gel(z,i) = gen_0;
     827    13408628 :   return z;
     828             : }
     829             : GEN
     830     1355666 : Rg_to_RgC(GEN x, long N)
     831             : {
     832     1355666 :   return (typ(x) == t_POL)? RgX_to_RgC(x,N): scalarcol_shallow(x, N);
     833             : }
     834             : 
     835             : /* vector of polynomials (in v) whose coefs are given by the columns of x */
     836             : GEN
     837      302131 : RgM_to_RgXV(GEN x, long v)
     838     1285114 : { pari_APPLY_type(t_VEC, RgV_to_RgX(gel(x,i), v)) }
     839             : GEN
     840        7202 : RgM_to_RgXV_reverse(GEN x, long v)
     841       28808 : { pari_APPLY_type(t_VEC, RgV_to_RgX_reverse(gel(x,i), v)) }
     842             : 
     843             : /* matrix whose entries are given by the coeffs of the polynomials in
     844             :  * vector v (considered as degree n-1 polynomials) */
     845             : GEN
     846      335749 : RgV_to_RgM(GEN x, long n)
     847     1686882 : { pari_APPLY_type(t_MAT, Rg_to_RgC(gel(x,i), n)) }
     848             : 
     849             : GEN
     850       78583 : RgXV_to_RgM(GEN x, long n)
     851      399634 : { pari_APPLY_type(t_MAT, RgX_to_RgC(gel(x,i), n)) }
     852             : 
     853             : /* polynomial (in v) of polynomials (in w) whose coeffs are given by the columns of x */
     854             : GEN
     855       23487 : RgM_to_RgXX(GEN x, long v,long w)
     856             : {
     857       23487 :   long j, lx = lg(x);
     858       23487 :   GEN y = cgetg(lx+1, t_POL);
     859       23487 :   y[1] = evalsigne(1) | evalvarn(v);
     860       23487 :   y++;
     861      130516 :   for (j=1; j<lx; j++) gel(y,j) = RgV_to_RgX(gel(x,j), w);
     862       23487 :   return normalizepol_lg(--y, lx+1);
     863             : }
     864             : 
     865             : /* matrix whose entries are given by the coeffs of the polynomial v in
     866             :  * two variables (considered as degree n-1 polynomials) */
     867             : GEN
     868         322 : RgXX_to_RgM(GEN v, long n)
     869             : {
     870         322 :   long j, N = lg(v)-1;
     871         322 :   GEN y = cgetg(N, t_MAT);
     872        1043 :   for (j=1; j<N; j++) gel(y,j) = Rg_to_RgC(gel(v,j+1), n);
     873         322 :   return y;
     874             : }
     875             : 
     876             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     877             : GEN
     878       32767 : RgXY_swapspec(GEN x, long n, long w, long nx)
     879             : {
     880       32767 :   long j, ly = n+3;
     881       32767 :   GEN y = cgetg(ly, t_POL);
     882       32767 :   y[1] = evalsigne(1);
     883      401015 :   for (j=2; j<ly; j++)
     884             :   {
     885             :     long k;
     886      368248 :     GEN a = cgetg(nx+2,t_POL);
     887      368248 :     a[1] = evalsigne(1) | evalvarn(w);
     888     2028272 :     for (k=0; k<nx; k++)
     889             :     {
     890     1660024 :       GEN xk = gel(x,k);
     891     1660024 :       if (typ(xk)==t_POL && varn(xk)==w)
     892     1565649 :         gel(a,k+2) = j<lg(xk)? gel(xk,j): gen_0;
     893             :       else
     894       94375 :         gel(a,k+2) = j==2 ? xk: gen_0;
     895             :     }
     896      368248 :     gel(y,j) = normalizepol_lg(a, nx+2);
     897             :   }
     898       32767 :   return normalizepol_lg(y,ly);
     899             : }
     900             : 
     901             : /* P(X,Y) --> P(Y,X), n is an upper bound for deg_Y(P) */
     902             : GEN
     903        2051 : RgXY_swap(GEN x, long n, long w)
     904             : {
     905        2051 :   GEN z = RgXY_swapspec(x+2, n, w, lgpol(x));
     906        2051 :   setvarn(z, varn(x)); return z;
     907             : }
     908             : 
     909             : long
     910        1387 : RgXY_degreex(GEN b)
     911             : {
     912        1387 :   long deg = 0, i;
     913        1387 :   if (!signe(b)) return -1;
     914       15592 :   for (i = 2; i < lg(b); ++i)
     915             :   {
     916       14205 :     GEN bi = gel(b, i);
     917       14205 :     if (typ(bi) == t_POL)
     918       13714 :       deg = maxss(deg, degpol(bi));
     919             :   }
     920        1387 :   return deg;
     921             : }
     922             : 
     923             : GEN
     924       38738 : RgXY_derivx(GEN x) { pari_APPLY_pol(RgX_deriv(gel(x,i))); }
     925             : 
     926             : /* return (x % X^n). Shallow */
     927             : GEN
     928     7937817 : RgXn_red_shallow(GEN a, long n)
     929             : {
     930     7937817 :   long i, L = n+2, l = lg(a);
     931             :   GEN  b;
     932     7937817 :   if (L >= l) return a; /* deg(x) < n */
     933     5748422 :   b = cgetg(L, t_POL); b[1] = a[1];
     934    37223148 :   for (i=2; i<L; i++) gel(b,i) = gel(a,i);
     935     5748423 :   return normalizepol_lg(b,L);
     936             : }
     937             : 
     938             : GEN
     939         483 : RgXnV_red_shallow(GEN x, long n)
     940        2268 : { pari_APPLY_type(t_VEC, RgXn_red_shallow(gel(x,i), n)) }
     941             : 
     942             : /* return (x * X^n). Shallow */
     943             : GEN
     944   175620096 : RgX_shift_shallow(GEN a, long n)
     945             : {
     946   175620096 :   long i, l = lg(a);
     947             :   GEN  b;
     948   175620096 :   if (l == 2 || !n) return a;
     949   110567113 :   l += n;
     950   110567113 :   if (n < 0)
     951             :   {
     952    55589122 :     if (l <= 2) return pol_0(varn(a));
     953    54117454 :     b = cgetg(l, t_POL); b[1] = a[1];
     954    54118221 :     a -= n;
     955   166557219 :     for (i=2; i<l; i++) gel(b,i) = gel(a,i);
     956             :   } else {
     957    54977991 :     b = cgetg(l, t_POL); b[1] = a[1];
     958    54982140 :     a -= n; n += 2;
     959   119498320 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     960   212046594 :     for (   ; i<l; i++) gel(b,i) = gel(a,i);
     961             :   }
     962   109100361 :   return b;
     963             : }
     964             : /* return (x * X^n). */
     965             : GEN
     966     1582243 : RgX_shift(GEN a, long n)
     967             : {
     968     1582243 :   long i, l = lg(a);
     969             :   GEN  b;
     970     1582243 :   if (l == 2 || !n) return RgX_copy(a);
     971     1581291 :   l += n;
     972     1581291 :   if (n < 0)
     973             :   {
     974        1442 :     if (l <= 2) return pol_0(varn(a));
     975        1372 :     b = cgetg(l, t_POL); b[1] = a[1];
     976        1372 :     a -= n;
     977        9534 :     for (i=2; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     978             :   } else {
     979     1579849 :     b = cgetg(l, t_POL); b[1] = a[1];
     980     1579849 :     a -= n; n += 2;
     981     3994994 :     for (i=2; i<n; i++) gel(b,i) = gen_0;
     982     4356380 :     for (   ; i<l; i++) gel(b,i) = gcopy(gel(a,i));
     983             :   }
     984     1581221 :   return b;
     985             : }
     986             : 
     987             : GEN
     988      317037 : RgX_rotate_shallow(GEN P, long k, long p)
     989             : {
     990      317037 :   long i, l = lgpol(P);
     991             :   GEN r;
     992      317037 :   if (signe(P)==0)
     993        1365 :     return pol_0(varn(P));
     994      315672 :   r = cgetg(p+2,t_POL); r[1] = P[1];
     995     2100644 :   for(i=0; i<p; i++)
     996             :   {
     997     1784972 :     long s = 2+(i+k)%p;
     998     1784972 :     gel(r,s) = i<l? gel(P,2+i): gen_0;
     999             :   }
    1000      315672 :   return RgX_renormalize(r);
    1001             : }
    1002             : 
    1003             : GEN
    1004     3000926 : RgX_mulXn(GEN x, long d)
    1005             : {
    1006             :   pari_sp av;
    1007             :   GEN z;
    1008             :   long v;
    1009     3000926 :   if (d >= 0) return RgX_shift(x, d);
    1010     1471950 :   d = -d;
    1011     1471950 :   v = RgX_val(x);
    1012     1471950 :   if (v >= d) return RgX_shift(x, -d);
    1013     1471936 :   av = avma;
    1014     1471936 :   z = gred_rfrac_simple(RgX_shift_shallow(x, -v), pol_xn(d - v, varn(x)));
    1015     1471936 :   return gerepileupto(av, z);
    1016             : }
    1017             : 
    1018             : long
    1019         588 : RgXV_maxdegree(GEN x)
    1020             : {
    1021         588 :   long d = -1, i, l = lg(x);
    1022        4494 :   for (i = 1; i < l; i++)
    1023        3906 :     d = maxss(d, degpol(gel(x,i)));
    1024         588 :   return d;
    1025             : }
    1026             : 
    1027             : long
    1028     3667177 : RgX_val(GEN x)
    1029             : {
    1030     3667177 :   long i, lx = lg(x);
    1031     3667177 :   if (lx == 2) return LONG_MAX;
    1032     4618932 :   for (i = 2; i < lx; i++)
    1033     4618876 :     if (!isexactzero(gel(x,i))) break;
    1034     3656250 :   if (i == lx) return LONG_MAX;/* possible with nonrational zeros */
    1035     3656194 :   return i - 2;
    1036             : }
    1037             : long
    1038    82477896 : RgX_valrem(GEN x, GEN *Z)
    1039             : {
    1040    82477896 :   long v, i, lx = lg(x);
    1041    82477896 :   if (lx == 2) { *Z = pol_0(varn(x)); return LONG_MAX; }
    1042   128451142 :   for (i = 2; i < lx; i++)
    1043   128452256 :     if (!isexactzero(gel(x,i))) break;
    1044             :   /* possible with nonrational zeros */
    1045    82477623 :   if (i == lx)
    1046             :   {
    1047          14 :     *Z = scalarpol_shallow(Rg_get_0(x), varn(x));
    1048          14 :     return LONG_MAX;
    1049             :   }
    1050    82477609 :   v = i - 2;
    1051    82477609 :   *Z = RgX_shift_shallow(x, -v);
    1052    82488994 :   return v;
    1053             : }
    1054             : long
    1055      852218 : RgX_valrem_inexact(GEN x, GEN *Z)
    1056             : {
    1057             :   long v;
    1058      852218 :   if (!signe(x)) { if (Z) *Z = pol_0(varn(x)); return LONG_MAX; }
    1059      869198 :   for (v = 0;; v++)
    1060      869198 :     if (!gequal0(gel(x,2+v))) break;
    1061      852211 :   if (Z) *Z = RgX_shift_shallow(x, -v);
    1062      852211 :   return v;
    1063             : }
    1064             : 
    1065             : GEN
    1066       67683 : RgXQC_red(GEN x, GEN T)
    1067      412979 : { pari_APPLY_type(t_COL, grem(gel(x,i), T)) }
    1068             : 
    1069             : GEN
    1070        1211 : RgXQV_red(GEN x, GEN T)
    1071       27832 : { pari_APPLY_type(t_VEC, grem(gel(x,i), T)) }
    1072             : 
    1073             : GEN
    1074       12999 : RgXQM_red(GEN x, GEN T)
    1075       80682 : { pari_APPLY_same(RgXQC_red(gel(x,i), T)) }
    1076             : 
    1077             : GEN
    1078         322 : RgXQM_mul(GEN P, GEN Q, GEN T)
    1079             : {
    1080         322 :   return RgXQM_red(RgM_mul(P, Q), T);
    1081             : }
    1082             : 
    1083             : GEN
    1084      507647 : RgXQX_red(GEN P, GEN T)
    1085             : {
    1086      507647 :   long i, l = lg(P);
    1087      507647 :   GEN Q = cgetg(l, t_POL);
    1088      507647 :   Q[1] = P[1];
    1089     2626850 :   for (i=2; i<l; i++) gel(Q,i) = grem(gel(P,i), T);
    1090      507647 :   return normalizepol_lg(Q, l);
    1091             : }
    1092             : 
    1093             : GEN
    1094      825589 : RgX_deriv(GEN x)
    1095             : {
    1096      825589 :   long i,lx = lg(x)-1;
    1097             :   GEN y;
    1098             : 
    1099      825589 :   if (lx<3) return pol_0(varn(x));
    1100      822502 :   y = cgetg(lx,t_POL); gel(y,2) = gcopy(gel(x,3));
    1101     3640771 :   for (i=3; i<lx ; i++) gel(y,i) = gmulsg(i-1,gel(x,i+1));
    1102      822492 :   y[1] = x[1]; return normalizepol_lg(y,i);
    1103             : }
    1104             : 
    1105             : GEN
    1106     2584304 : RgX_recipspec_shallow(GEN x, long l, long n)
    1107             : {
    1108             :   long i;
    1109     2584304 :   GEN z = cgetg(n+2,t_POL);
    1110     2584311 :   z[1] = 0; z += 2;
    1111   144282596 :   for(i=0; i<l; i++) gel(z,n-i-1) = gel(x,i);
    1112     2834082 :   for(   ; i<n; i++) gel(z, n-i-1) = gen_0;
    1113     2584311 :   return normalizepol_lg(z-2,n+2);
    1114             : }
    1115             : 
    1116             : GEN
    1117      643633 : RgXn_recip_shallow(GEN P, long n)
    1118             : {
    1119      643633 :   GEN Q = RgX_recipspec_shallow(P+2, lgpol(P), n);
    1120      643641 :   setvarn(Q, varn(P));
    1121      643641 :   return Q;
    1122             : }
    1123             : 
    1124             : /* return coefficients s.t x = x_0 X^n + ... + x_n */
    1125             : GEN
    1126       31367 : RgX_recip(GEN x)
    1127             : {
    1128             :   long lx, i, j;
    1129       31367 :   GEN y = cgetg_copy(x, &lx);
    1130      274484 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gcopy(gel(x,j));
    1131       31367 :   return normalizepol_lg(y,lx);
    1132             : }
    1133             : /* shallow version */
    1134             : GEN
    1135       59388 : RgX_recip_shallow(GEN x)
    1136             : {
    1137             :   long lx, i, j;
    1138       59388 :   GEN y = cgetg_copy(x, &lx);
    1139      356118 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
    1140       59388 :   return normalizepol_lg(y,lx);
    1141             : }
    1142             : 
    1143             : GEN
    1144     5120233 : RgX_recip_i(GEN x)
    1145             : {
    1146             :   long lx, i, j;
    1147     5120233 :   GEN y = cgetg_copy(x, &lx);
    1148    25977514 :   y[1] = x[1]; for (i=2,j=lx-1; i<lx; i++,j--) gel(y,i) = gel(x,j);
    1149     5120240 :   return y;
    1150             : }
    1151             : /*******************************************************************/
    1152             : /*                                                                 */
    1153             : /*                      ADDITION / SUBTRACTION                     */
    1154             : /*                                                                 */
    1155             : /*******************************************************************/
    1156             : /* same variable */
    1157             : GEN
    1158   145062912 : RgX_add(GEN x, GEN y)
    1159             : {
    1160   145062912 :   long i, lx = lg(x), ly = lg(y);
    1161             :   GEN z;
    1162   145062912 :   if (ly <= lx) {
    1163   130827333 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1164   509071349 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1165   188651617 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1166   130802152 :     z = normalizepol_lg(z, lx);
    1167             :   } else {
    1168    14235579 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1169    54267195 :     for (i=2; i < lx; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
    1170    41660010 :     for (   ; i < ly; i++) gel(z,i) = gcopy(gel(y,i));
    1171    14236848 :     z = normalizepol_lg(z, ly);
    1172             :   }
    1173   145052987 :   return z;
    1174             : }
    1175             : GEN
    1176    71038154 : RgX_sub(GEN x, GEN y)
    1177             : {
    1178    71038154 :   long i, lx = lg(x), ly = lg(y);
    1179             :   GEN z;
    1180    71038154 :   if (ly <= lx) {
    1181    37178271 :     z = cgetg(lx,t_POL); z[1] = x[1];
    1182   185049236 :     for (i=2; i < ly; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1183    63749105 :     for (   ; i < lx; i++) gel(z,i) = gcopy(gel(x,i));
    1184    37154772 :     z = normalizepol_lg(z, lx);
    1185             :   } else {
    1186    33859883 :     z = cgetg(ly,t_POL); z[1] = y[1];
    1187   125835794 :     for (i=2; i < lx; i++) gel(z,i) = gsub(gel(x,i),gel(y,i));
    1188    74990413 :     for (   ; i < ly; i++) gel(z,i) = gneg(gel(y,i));
    1189    33828333 :     z = normalizepol_lg(z, ly);
    1190             :   }
    1191    71036371 :   return z;
    1192             : }
    1193             : GEN
    1194     8837971 : RgX_neg(GEN x)
    1195    53209269 : { pari_APPLY_pol_normalized(gneg(gel(x,i))); }
    1196             : 
    1197             : GEN
    1198    44558189 : RgX_Rg_add(GEN y, GEN x)
    1199             : {
    1200             :   GEN z;
    1201    44558189 :   long lz = lg(y), i;
    1202    44558189 :   if (lz == 2) return scalarpol(x,varn(y));
    1203    22669975 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1204    22669948 :   gel(z,2) = gadd(gel(y,2),x);
    1205    78843802 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1206             :   /* probably useless unless lz = 3, but cannot be skipped if y is
    1207             :    * an inexact 0 */
    1208    22669869 :   return normalizepol_lg(z,lz);
    1209             : }
    1210             : GEN
    1211       65227 : RgX_Rg_add_shallow(GEN y, GEN x)
    1212             : {
    1213             :   GEN z;
    1214       65227 :   long lz = lg(y), i;
    1215       65227 :   if (lz == 2) return scalarpol(x,varn(y));
    1216       65227 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1217       65227 :   gel(z,2) = gadd(gel(y,2),x);
    1218      130614 :   for(i=3; i<lz; i++) gel(z,i) = gel(y,i);
    1219       65226 :   return normalizepol_lg(z,lz);
    1220             : }
    1221             : GEN
    1222     1223134 : RgX_Rg_sub(GEN y, GEN x)
    1223             : {
    1224             :   GEN z;
    1225     1223134 :   long lz = lg(y), i;
    1226     1223134 :   if (lz == 2)
    1227             :   { /* scalarpol(gneg(x),varn(y)) optimized */
    1228      922321 :     long v = varn(y);
    1229      922321 :     if (isrationalzero(x)) return pol_0(v);
    1230        2235 :     z = cgetg(3,t_POL);
    1231        2235 :     z[1] = gequal0(x)? evalvarn(v)
    1232        2235 :                    : evalvarn(v) | evalsigne(1);
    1233        2235 :     gel(z,2) = gneg(x); return z;
    1234             :   }
    1235      300813 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1236      300813 :   gel(z,2) = gsub(gel(y,2),x);
    1237      689236 :   for(i=3; i<lz; i++) gel(z,i) = gcopy(gel(y,i));
    1238      300813 :   return normalizepol_lg(z,lz);
    1239             : }
    1240             : GEN
    1241     4952191 : Rg_RgX_sub(GEN x, GEN y)
    1242             : {
    1243             :   GEN z;
    1244     4952191 :   long lz = lg(y), i;
    1245     4952191 :   if (lz == 2) return scalarpol(x,varn(y));
    1246     4882196 :   z = cgetg(lz,t_POL); z[1] = y[1];
    1247     4882193 :   gel(z,2) = gsub(x, gel(y,2));
    1248     7329262 :   for(i=3; i<lz; i++) gel(z,i) = gneg(gel(y,i));
    1249     4882139 :   return normalizepol_lg(z,lz);
    1250             : }
    1251             : /*******************************************************************/
    1252             : /*                                                                 */
    1253             : /*                  KARATSUBA MULTIPLICATION                       */
    1254             : /*                                                                 */
    1255             : /*******************************************************************/
    1256             : #if 0
    1257             : /* to debug Karatsuba-like routines */
    1258             : GEN
    1259             : zx_debug_spec(GEN x, long nx)
    1260             : {
    1261             :   GEN z = cgetg(nx+2,t_POL);
    1262             :   long i;
    1263             :   for (i=0; i<nx; i++) gel(z,i+2) = stoi(x[i]);
    1264             :   z[1] = evalsigne(1); return z;
    1265             : }
    1266             : 
    1267             : GEN
    1268             : RgX_debug_spec(GEN x, long nx)
    1269             : {
    1270             :   GEN z = cgetg(nx+2,t_POL);
    1271             :   long i;
    1272             :   for (i=0; i<nx; i++) z[i+2] = x[i];
    1273             :   z[1] = evalsigne(1); return z;
    1274             : }
    1275             : #endif
    1276             : 
    1277             : /* generic multiplication */
    1278             : GEN
    1279     8912166 : RgX_addspec_shallow(GEN x, GEN y, long nx, long ny)
    1280             : {
    1281             :   GEN z, t;
    1282             :   long i;
    1283     8912166 :   if (nx == ny) {
    1284     1535561 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1285     4952612 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1286     1535553 :     return normalizepol_lg(z, nx+2);
    1287             :   }
    1288     7376605 :   if (ny < nx) {
    1289     7189906 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1290    26470345 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1291    17182156 :     for (   ; i < nx; i++) gel(t,i) = gel(x,i);
    1292     7189441 :     return normalizepol_lg(z, nx+2);
    1293             :   } else {
    1294      186699 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1295     3678005 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1296      451698 :     for (   ; i < ny; i++) gel(t,i) = gel(y,i);
    1297      186715 :     return normalizepol_lg(z, ny+2);
    1298             :   }
    1299             : }
    1300             : GEN
    1301      223426 : RgX_addspec(GEN x, GEN y, long nx, long ny)
    1302             : {
    1303             :   GEN z, t;
    1304             :   long i;
    1305      223426 :   if (nx == ny) {
    1306       12824 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1307     2185778 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1308       12824 :     return normalizepol_lg(z, nx+2);
    1309             :   }
    1310      210602 :   if (ny < nx) {
    1311      208824 :     z = cgetg(nx+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1312     3729526 :     for (i=0; i < ny; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1313     2382089 :     for (   ; i < nx; i++) gel(t,i) = gcopy(gel(x,i));
    1314      208824 :     return normalizepol_lg(z, nx+2);
    1315             :   } else {
    1316        1778 :     z = cgetg(ny+2,t_POL); z[1] = evalsigne(1)|evalvarn(0); t = z+2;
    1317      330904 :     for (i=0; i < nx; i++) gel(t,i) = gadd(gel(x,i),gel(y,i));
    1318       12222 :     for (   ; i < ny; i++) gel(t,i) = gcopy(gel(y,i));
    1319        1778 :     return normalizepol_lg(z, ny+2);
    1320             :   }
    1321             : }
    1322             : 
    1323             : /* Return the vector of coefficients of x, where we replace rational 0s by NULL
    1324             :  * [ to speed up basic operation s += x[i]*y[j] ]. We create a proper
    1325             :  * t_VECSMALL, to hold this, which can be left on stack: gerepile
    1326             :  * will not crash on it. The returned vector itself is not a proper GEN,
    1327             :  * we access the coefficients as x[i], i = 0..deg(x) */
    1328             : static GEN
    1329    77672523 : RgXspec_kill0(GEN x, long lx)
    1330             : {
    1331    77672523 :   GEN z = cgetg(lx+1, t_VECSMALL) + 1; /* inhibit gerepile-wise */
    1332             :   long i;
    1333   238712177 :   for (i=0; i <lx; i++)
    1334             :   {
    1335   161039656 :     GEN c = gel(x,i);
    1336   161039656 :     z[i] = (long)(isrationalzero(c)? NULL: c);
    1337             :   }
    1338    77672521 :   return z;
    1339             : }
    1340             : 
    1341             : INLINE GEN
    1342   109230118 : RgX_mulspec_basecase_limb(GEN x, GEN y, long a, long b)
    1343             : {
    1344   109230118 :   pari_sp av = avma;
    1345   109230118 :   GEN s = NULL;
    1346             :   long i;
    1347             : 
    1348   388409068 :   for (i=a; i<b; i++)
    1349   279188414 :     if (gel(y,i) && gel(x,-i))
    1350             :     {
    1351   200386466 :       GEN t = gmul(gel(y,i), gel(x,-i));
    1352   200381889 :       s = s? gadd(s, t): t;
    1353             :     }
    1354   109220654 :   return s? gerepileupto(av, s): gen_0;
    1355             : }
    1356             : 
    1357             : /* assume nx >= ny > 0, return x * y * t^v */
    1358             : static GEN
    1359    32189590 : RgX_mulspec_basecase(GEN x, GEN y, long nx, long ny, long v)
    1360             : {
    1361             :   long i, lz, nz;
    1362             :   GEN z;
    1363             : 
    1364    32189590 :   x = RgXspec_kill0(x,nx);
    1365    32189549 :   y = RgXspec_kill0(y,ny);
    1366    32189567 :   lz = nx + ny + 1; nz = lz-2;
    1367    32189567 :   lz += v;
    1368    32189567 :   z = cgetg(lz, t_POL) + 2; /* x:y:z [i] = term of degree i */
    1369    69227113 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1370    79223922 :   for (i=0; i<ny; i++)gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0, i+1);
    1371    54700007 :   for (  ; i<nx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ny);
    1372    47034500 :   for (  ; i<nz; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-nx+1,ny);
    1373    32188997 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1374             : }
    1375             : 
    1376             : /* return (x * X^d) + y. Assume d > 0 */
    1377             : GEN
    1378    10068500 : RgX_addmulXn_shallow(GEN x0, GEN y0, long d)
    1379             : {
    1380             :   GEN x, y, xd, yd, zd;
    1381             :   long a, lz, nx, ny;
    1382             : 
    1383    10068500 :   if (!signe(x0)) return y0;
    1384     9433906 :   ny = lgpol(y0);
    1385     9433900 :   nx = lgpol(x0);
    1386     9434008 :   zd = (GEN)avma;
    1387     9434008 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1388     9434008 :   if (a <= 0)
    1389             :   {
    1390     1479995 :     lz = nx+d+2;
    1391     1479995 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1392     3398604 :     while (xd > x) gel(--zd,0) = gel(--xd,0);
    1393     1479994 :     x = zd + a;
    1394     1497702 :     while (zd > x) gel(--zd,0) = gen_0;
    1395             :   }
    1396             :   else
    1397             :   {
    1398     7954013 :     xd = new_chunk(d); yd = y+d;
    1399     7954020 :     x = RgX_addspec_shallow(x,yd, nx,a);
    1400     7953952 :     lz = (a>nx)? ny+2: lg(x)+d;
    1401    39643700 :     x += 2; while (xd > x) *--zd = *--xd;
    1402             :   }
    1403    22984412 :   while (yd > y) *--zd = *--yd;
    1404     9433946 :   *--zd = x0[1];
    1405     9433946 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1406             : }
    1407             : GEN
    1408      515999 : RgX_addmulXn(GEN x0, GEN y0, long d)
    1409             : {
    1410             :   GEN x, y, xd, yd, zd;
    1411             :   long a, lz, nx, ny;
    1412             : 
    1413      515999 :   if (!signe(x0)) return RgX_copy(y0);
    1414      515215 :   nx = lgpol(x0);
    1415      515215 :   ny = lgpol(y0);
    1416      515215 :   zd = (GEN)avma;
    1417      515215 :   x = x0 + 2; y = y0 + 2; a = ny-d;
    1418      515215 :   if (a <= 0)
    1419             :   {
    1420      291789 :     lz = nx+d+2;
    1421      291789 :     (void)new_chunk(lz); xd = x+nx; yd = y+ny;
    1422     4294362 :     while (xd > x) gel(--zd,0) = gcopy(gel(--xd,0));
    1423      291789 :     x = zd + a;
    1424      757955 :     while (zd > x) gel(--zd,0) = gen_0;
    1425             :   }
    1426             :   else
    1427             :   {
    1428      223426 :     xd = new_chunk(d); yd = y+d;
    1429      223426 :     x = RgX_addspec(x,yd, nx,a);
    1430      223426 :     lz = (a>nx)? ny+2: lg(x)+d;
    1431     8429917 :     x += 2; while (xd > x) *--zd = *--xd;
    1432             :   }
    1433     2622808 :   while (yd > y) gel(--zd,0) = gcopy(gel(--yd,0));
    1434      515215 :   *--zd = x0[1];
    1435      515215 :   *--zd = evaltyp(t_POL) | evallg(lz); return zd;
    1436             : }
    1437             : 
    1438             : /* return x * y mod t^n */
    1439             : static GEN
    1440     6627959 : RgXn_mul_basecase(GEN x, GEN y, long n)
    1441             : {
    1442     6627959 :   long i, lz = n+2, lx = lgpol(x), ly = lgpol(y);
    1443             :   GEN z;
    1444     6627959 :   if (lx < 0) return pol_0(varn(x));
    1445     6627959 :   if (ly < 0) return pol_0(varn(x));
    1446     6627959 :   z = cgetg(lz, t_POL) + 2;
    1447     6627959 :   x+=2; if (lx > n) lx = n;
    1448     6627959 :   y+=2; if (ly > n) ly = n;
    1449     6627959 :   z[-1] = x[-1];
    1450     6627959 :   if (ly > lx) { swap(x,y); lswap(lx,ly); }
    1451     6627959 :   x = RgXspec_kill0(x, lx);
    1452     6627959 :   y = RgXspec_kill0(y, ly);
    1453             :   /* x:y:z [i] = term of degree i */
    1454    26223885 :   for (i=0;i<ly; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,i+1);
    1455    11826184 :   for (  ; i<lx; i++) gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, 0,ly);
    1456     6674031 :   for (  ; i<n; i++)  gel(z,i) = RgX_mulspec_basecase_limb(x+i,y, i-lx+1,ly);
    1457     6627959 :   return normalizepol_lg(z - 2, lz);
    1458             : }
    1459             : /* Mulders / Karatsuba product f*g mod t^n (Hanrot-Zimmermann variant) */
    1460             : static GEN
    1461     9413940 : RgXn_mul2(GEN f, GEN g, long n)
    1462             : {
    1463     9413940 :   pari_sp av = avma;
    1464             :   GEN fe,fo, ge,go, l,h,m;
    1465             :   long n0, n1;
    1466     9413940 :   if (degpol(f) + degpol(g) < n) return RgX_mul(f,g);
    1467     6663449 :   if (n < 80) return RgXn_mul_basecase(f,g,n);
    1468       35490 :   n0 = n>>1; n1 = n-n0;
    1469       35490 :   RgX_even_odd(f, &fe, &fo);
    1470       35490 :   RgX_even_odd(g, &ge, &go);
    1471       35490 :   l = RgXn_mul2(fe,ge,n1);
    1472       35490 :   h = RgXn_mul2(fo,go,n0);
    1473       35490 :   m = RgX_sub(RgXn_mul2(RgX_add(fe,fo),RgX_add(ge,go),n0), RgX_add(l,h));
    1474             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1475             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1476       35490 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1477             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1478       35490 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1479       35490 :   m = RgX_inflate(m,2);
    1480             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1481       35490 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1482       35490 :   h = RgX_inflate(h,2);
    1483       35490 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1484       35490 :   return gerepileupto(av, h);
    1485             : }
    1486             : /* (f*g) \/ x^n */
    1487             : static GEN
    1488     1589258 : RgX_mulhigh_i2(GEN f, GEN g, long n)
    1489             : {
    1490     1589258 :   long d = degpol(f)+degpol(g) + 1 - n;
    1491             :   GEN h;
    1492     1589258 :   if (d <= 2) return RgX_shift_shallow(RgX_mul(f,g), -n);
    1493       29654 :   h = RgX_recip_i(RgXn_mul2(RgX_recip_i(f),
    1494             :                                   RgX_recip_i(g), d));
    1495       29654 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1496             : }
    1497             : 
    1498             : /* (f*g) \/ x^n */
    1499             : static GEN
    1500           0 : RgX_sqrhigh_i2(GEN f, long n)
    1501             : {
    1502           0 :   long d = 2*degpol(f)+ 1 - n;
    1503             :   GEN h;
    1504           0 :   if (d <= 2) return RgX_shift_shallow(RgX_sqr(f), -n);
    1505           0 :   h = RgX_recip_i(RgXn_sqr(RgX_recip_i(f), d));
    1506           0 :   return RgX_shift_shallow(h, d-1-degpol(h)); /* possibly (fg)(0) = 0 */
    1507             : }
    1508             : 
    1509             : /* fast product (Karatsuba) of polynomials a,b. These are not real GENs, a+2,
    1510             :  * b+2 were sent instead. na, nb = number of terms of a, b.
    1511             :  * Only c, c0, c1, c2 are genuine GEN.
    1512             :  */
    1513             : GEN
    1514    38252588 : RgX_mulspec(GEN a, GEN b, long na, long nb)
    1515             : {
    1516             :   GEN a0, c, c0;
    1517    38252588 :   long n0, n0a, i, v = 0;
    1518             :   pari_sp av;
    1519             : 
    1520    63516548 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v++; }
    1521    58136920 :   while (nb && isrationalzero(gel(b,0))) { b++; nb--; v++; }
    1522    38252571 :   if (na < nb) swapspec(a,b, na,nb);
    1523    38252571 :   if (!nb) return pol_0(0);
    1524             : 
    1525    32670047 :   if (nb < RgX_MUL_LIMIT) return RgX_mulspec_basecase(a,b,na,nb, v);
    1526      480464 :   RgX_shift_inplace_init(v);
    1527      480470 :   i = (na>>1); n0 = na-i; na = i;
    1528      480470 :   av = avma; a0 = a+n0; n0a = n0;
    1529     1335986 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1530             : 
    1531      480470 :   if (nb > n0)
    1532             :   {
    1533             :     GEN b0,c1,c2;
    1534             :     long n0b;
    1535             : 
    1536      479071 :     nb -= n0; b0 = b+n0; n0b = n0;
    1537     1447050 :     while (n0b && isrationalzero(gel(b,n0b-1))) n0b--;
    1538      479071 :     c = RgX_mulspec(a,b,n0a,n0b);
    1539      479071 :     c0 = RgX_mulspec(a0,b0, na,nb);
    1540             : 
    1541      479071 :     c2 = RgX_addspec_shallow(a0,a, na,n0a);
    1542      479071 :     c1 = RgX_addspec_shallow(b0,b, nb,n0b);
    1543             : 
    1544      479071 :     c1 = RgX_mulspec(c1+2,c2+2, lgpol(c1),lgpol(c2));
    1545      479071 :     c2 = RgX_sub(c1, RgX_add(c0,c));
    1546      479071 :     c0 = RgX_addmulXn_shallow(c0, c2, n0);
    1547             :   }
    1548             :   else
    1549             :   {
    1550        1399 :     c = RgX_mulspec(a,b,n0a,nb);
    1551        1399 :     c0 = RgX_mulspec(a0,b,na,nb);
    1552             :   }
    1553      480470 :   c0 = RgX_addmulXn(c0,c,n0);
    1554      480470 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1555             : }
    1556             : 
    1557             : INLINE GEN
    1558       98675 : RgX_sqrspec_basecase_limb(GEN x, long a, long i)
    1559             : {
    1560       98675 :   pari_sp av = avma;
    1561       98675 :   GEN s = NULL;
    1562       98675 :   long j, l = (i+1)>>1;
    1563      201607 :   for (j=a; j<l; j++)
    1564             :   {
    1565      102932 :     GEN xj = gel(x,j), xx = gel(x,i-j);
    1566      102932 :     if (xj && xx)
    1567             :     {
    1568       93846 :       GEN t = gmul(xj, xx);
    1569       93846 :       s = s? gadd(s, t): t;
    1570             :     }
    1571             :   }
    1572       98675 :   if (s) s = gshift(s,1);
    1573       98675 :   if ((i&1) == 0)
    1574             :   {
    1575       68152 :     GEN t = gel(x, i>>1);
    1576       68152 :     if (t) {
    1577       65107 :       t = gsqr(t);
    1578       65107 :       s = s? gadd(s, t): t;
    1579             :     }
    1580             :   }
    1581       98675 :   return s? gerepileupto(av,s): gen_0;
    1582             : }
    1583             : static GEN
    1584       37629 : RgX_sqrspec_basecase(GEN x, long nx, long v)
    1585             : {
    1586             :   long i, lz, nz;
    1587             :   GEN z;
    1588             : 
    1589       37629 :   if (!nx) return pol_0(0);
    1590       37629 :   x = RgXspec_kill0(x,nx);
    1591       37629 :   lz = (nx << 1) + 1, nz = lz-2;
    1592       37629 :   lz += v;
    1593       37629 :   z = cgetg(lz,t_POL) + 2;
    1594       93727 :   for (i=0; i<v; i++) gel(z++, 0) = gen_0;
    1595      105781 :   for (i=0; i<nx; i++)gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1596       68152 :   for (  ; i<nz; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, i-nx+1, i);
    1597       37629 :   z -= v+2; z[1] = 0; return normalizepol_lg(z, lz);
    1598             : }
    1599             : /* return x^2 mod t^n */
    1600             : static GEN
    1601           0 : RgXn_sqr_basecase(GEN x, long n)
    1602             : {
    1603           0 :   long i, lz = n+2, lx = lgpol(x);
    1604             :   GEN z;
    1605           0 :   if (lx < 0) return pol_0(varn(x));
    1606           0 :   z = cgetg(lz, t_POL);
    1607           0 :   z[1] = x[1];
    1608           0 :   x+=2; if (lx > n) lx = n;
    1609           0 :   x = RgXspec_kill0(x,lx);
    1610           0 :   z+=2;/* x:z [i] = term of degree i */
    1611           0 :   for (i=0;i<lx; i++) gel(z,i) = RgX_sqrspec_basecase_limb(x, 0, i);
    1612           0 :   for (  ; i<n; i++)  gel(z,i) = RgX_sqrspec_basecase_limb(x, i-lx+1, i);
    1613           0 :   z -= 2; return normalizepol_lg(z, lz);
    1614             : }
    1615             : /* Mulders / Karatsuba product f^2 mod t^n (Hanrot-Zimmermann variant) */
    1616             : static GEN
    1617         315 : RgXn_sqr2(GEN f, long n)
    1618             : {
    1619         315 :   pari_sp av = avma;
    1620             :   GEN fe,fo, l,h,m;
    1621             :   long n0, n1;
    1622         315 :   if (2*degpol(f) < n) return RgX_sqr_i(f);
    1623           0 :   if (n < 80) return RgXn_sqr_basecase(f,n);
    1624           0 :   n0 = n>>1; n1 = n-n0;
    1625           0 :   RgX_even_odd(f, &fe, &fo);
    1626           0 :   l = RgXn_sqr(fe,n1);
    1627           0 :   h = RgXn_sqr(fo,n0);
    1628           0 :   m = RgX_sub(RgXn_sqr(RgX_add(fe,fo),n0), RgX_add(l,h));
    1629             :   /* n1-1 <= n0 <= n1, deg l,m <= n1-1, deg h <= n0-1
    1630             :    * result is t^2 h(t^2) + t m(t^2) + l(t^2) */
    1631           0 :   l = RgX_inflate(l,2); /* deg l <= 2n1 - 2 <= n-1 */
    1632             :   /* deg(t m(t^2)) <= 2n1 - 1 <= n, truncate to < n */
    1633           0 :   if (2*degpol(m)+1 == n) m = normalizepol_lg(m, lg(m)-1);
    1634           0 :   m = RgX_inflate(m,2);
    1635             :   /* deg(t^2 h(t^2)) <= 2n0 <= n, truncate to < n */
    1636           0 :   if (2*degpol(h)+2 == n) h = normalizepol_lg(h, lg(h)-1);
    1637           0 :   h = RgX_inflate(h,2);
    1638           0 :   h = RgX_addmulXn(RgX_addmulXn_shallow(h,m,1), l,1);
    1639           0 :   return gerepileupto(av, h);
    1640             : }
    1641             : GEN
    1642       37668 : RgX_sqrspec(GEN a, long na)
    1643             : {
    1644             :   GEN a0, c, c0, c1;
    1645       37668 :   long n0, n0a, i, v = 0;
    1646             :   pari_sp av;
    1647             : 
    1648       65717 :   while (na && isrationalzero(gel(a,0))) { a++; na--; v += 2; }
    1649       37668 :   if (na<RgX_SQR_LIMIT) return RgX_sqrspec_basecase(a, na, v);
    1650          39 :   RgX_shift_inplace_init(v);
    1651          39 :   i = (na>>1); n0 = na-i; na = i;
    1652          39 :   av = avma; a0 = a+n0; n0a = n0;
    1653          39 :   while (n0a && isrationalzero(gel(a,n0a-1))) n0a--;
    1654             : 
    1655          39 :   c = RgX_sqrspec(a,n0a);
    1656          39 :   c0 = RgX_sqrspec(a0,na);
    1657          39 :   c1 = gmul2n(RgX_mulspec(a0,a, na,n0a), 1);
    1658          39 :   c0 = RgX_addmulXn_shallow(c0,c1, n0);
    1659          39 :   c0 = RgX_addmulXn(c0,c,n0);
    1660          39 :   return RgX_shift_inplace(gerepileupto(av,c0), v);
    1661             : }
    1662             : 
    1663             : /* (X^a + A)(X^b + B) - X^(a+b), where deg A < a, deg B < b */
    1664             : GEN
    1665     1719154 : RgX_mul_normalized(GEN A, long a, GEN B, long b)
    1666             : {
    1667     1719154 :   GEN z = RgX_mul(A, B);
    1668     1719130 :   if (a < b)
    1669       10442 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(A, B, b-a), z, a);
    1670     1708688 :   else if (a > b)
    1671     1104887 :     z = RgX_addmulXn_shallow(RgX_addmulXn_shallow(B, A, a-b), z, b);
    1672             :   else
    1673      603801 :     z = RgX_addmulXn_shallow(RgX_add(A, B), z, a);
    1674     1719146 :   return z;
    1675             : }
    1676             : 
    1677             : GEN
    1678    36812554 : RgX_mul_i(GEN x, GEN y)
    1679             : {
    1680    36812554 :   GEN z = RgX_mulspec(x+2, y+2, lgpol(x), lgpol(y));
    1681    36812185 :   setvarn(z, varn(x)); return z;
    1682             : }
    1683             : 
    1684             : GEN
    1685       37590 : RgX_sqr_i(GEN x)
    1686             : {
    1687       37590 :   GEN z = RgX_sqrspec(x+2, lgpol(x));
    1688       37590 :   setvarn(z,varn(x)); return z;
    1689             : }
    1690             : 
    1691             : /*******************************************************************/
    1692             : /*                                                                 */
    1693             : /*                               DIVISION                          */
    1694             : /*                                                                 */
    1695             : /*******************************************************************/
    1696             : GEN
    1697     7167749 : RgX_Rg_divexact(GEN x, GEN y) {
    1698     7167749 :   long i, lx = lg(x);
    1699             :   GEN z;
    1700     7167749 :   if (lx == 2) return gcopy(x);
    1701     7116959 :   switch(typ(y))
    1702             :   {
    1703     6956489 :     case t_INT:
    1704     6956489 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1705     6580460 :       break;
    1706        5243 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1707             :   }
    1708     6735687 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1709    34394830 :   for (i=2; i<lx; i++) gel(z,i) = gdivexact(gel(x,i),y);
    1710     6734006 :   return z;
    1711             : }
    1712             : GEN
    1713    30284753 : RgX_Rg_div(GEN x, GEN y) {
    1714    30284753 :   long i, lx = lg(x);
    1715             :   GEN z;
    1716    30284753 :   if (lx == 2) return gcopy(x);
    1717    30012009 :   switch(typ(y))
    1718             :   {
    1719    20784758 :     case t_INT:
    1720    20784758 :       if (is_pm1(y)) return signe(y) < 0 ? RgX_neg(x): RgX_copy(x);
    1721     4145669 :       break;
    1722        6244 :     case t_INTMOD: case t_POLMOD: return RgX_Rg_mul(x, ginv(y));
    1723             :   }
    1724    13366676 :   z = cgetg(lx, t_POL); z[1] = x[1];
    1725    49874446 :   for (i=2; i<lx; i++) gel(z,i) = gdiv(gel(x,i),y);
    1726    13366223 :   return normalizepol_lg(z, lx);
    1727             : }
    1728             : GEN
    1729       39130 : RgX_normalize(GEN x)
    1730             : {
    1731       39130 :   GEN z, d = NULL;
    1732       39130 :   long i, n = lg(x)-1;
    1733       39130 :   for (i = n; i > 1; i--) { d = gel(x,i); if (!gequal0(d)) break; }
    1734       39130 :   if (i == 1) return pol_0(varn(x));
    1735       39130 :   if (i == n && isint1(d)) return x;
    1736       17612 :   n = i; z = cgetg(n+1, t_POL); z[1] = x[1];
    1737       31612 :   for (i=2; i<n; i++) gel(z,i) = gdiv(gel(x,i),d);
    1738       17612 :   gel(z,n) = Rg_get_1(d); return z;
    1739             : }
    1740             : GEN
    1741       10458 : RgX_divs(GEN x, long y) { pari_APPLY_pol(gdivgs(gel(x,i),y)); }
    1742             : GEN
    1743      251522 : RgX_div_by_X_x(GEN a, GEN x, GEN *r)
    1744             : {
    1745      251522 :   long l = lg(a), i;
    1746             :   GEN a0, z0, z;
    1747             : 
    1748      251522 :   if (l <= 3)
    1749             :   {
    1750           0 :     if (r) *r = l == 2? gen_0: gcopy(gel(a,2));
    1751           0 :     return pol_0(varn(a));
    1752             :   }
    1753      251522 :   z = cgetg(l-1, t_POL);
    1754      251526 :   z[1] = a[1];
    1755      251526 :   a0 = a + l-1;
    1756      251526 :   z0 = z + l-2; *z0 = *a0--;
    1757     2956644 :   for (i=l-3; i>1; i--) /* z[i] = a[i+1] + x*z[i+1] */
    1758             :   {
    1759     2705123 :     GEN t = gadd(gel(a0--,0), gmul(x, gel(z0--,0)));
    1760     2705118 :     gel(z0,0) = t;
    1761             :   }
    1762      251521 :   if (r) *r = gadd(gel(a0,0), gmul(x, gel(z0,0)));
    1763      251521 :   return z;
    1764             : }
    1765             : /* Polynomial division x / y:
    1766             :  *   if pr = ONLY_REM return remainder, otherwise return quotient
    1767             :  *   if pr = ONLY_DIVIDES return quotient if division is exact, else NULL
    1768             :  *   if pr != NULL set *pr to remainder, as the last object on stack */
    1769             : /* assume, typ(x) = typ(y) = t_POL, same variable */
    1770             : static GEN
    1771    23519194 : RgX_divrem_i(GEN x, GEN y, GEN *pr)
    1772             : {
    1773             :   pari_sp avy, av, av1;
    1774             :   long dx,dy,dz,i,j,sx,lr;
    1775             :   GEN z,p1,p2,rem,y_lead,mod,p;
    1776             :   GEN (*f)(GEN,GEN);
    1777             : 
    1778    23519194 :   if (!signe(y)) pari_err_INV("RgX_divrem",y);
    1779             : 
    1780    23519194 :   dy = degpol(y);
    1781    23519172 :   y_lead = gel(y,dy+2);
    1782    23519172 :   if (gequal0(y_lead)) /* normalize denominator if leading term is 0 */
    1783             :   {
    1784           0 :     pari_warn(warner,"normalizing a polynomial with 0 leading term");
    1785           0 :     for (dy--; dy>=0; dy--)
    1786             :     {
    1787           0 :       y_lead = gel(y,dy+2);
    1788           0 :       if (!gequal0(y_lead)) break;
    1789             :     }
    1790             :   }
    1791    23519165 :   if (!dy) /* y is constant */
    1792             :   {
    1793        6857 :     if (pr == ONLY_REM) return pol_0(varn(x));
    1794        6850 :     z = RgX_Rg_div(x, y_lead);
    1795        6850 :     if (pr == ONLY_DIVIDES) return z;
    1796        2636 :     if (pr) *pr = pol_0(varn(x));
    1797        2636 :     return z;
    1798             :   }
    1799    23512308 :   dx = degpol(x);
    1800    23512277 :   if (dx < dy)
    1801             :   {
    1802     3850101 :     if (pr == ONLY_REM) return RgX_copy(x);
    1803      348542 :     if (pr == ONLY_DIVIDES) return signe(x)? NULL: pol_0(varn(x));
    1804      348521 :     z = pol_0(varn(x));
    1805      348521 :     if (pr) *pr = RgX_copy(x);
    1806      348521 :     return z;
    1807             :   }
    1808             : 
    1809             :   /* x,y in R[X], y non constant */
    1810    19662176 :   av = avma;
    1811    19662176 :   p = NULL;
    1812    19662176 :   if (RgX_is_FpX(x, &p) && RgX_is_FpX(y, &p) && p)
    1813             :   {
    1814      130403 :     z = FpX_divrem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p, pr);
    1815      130403 :     if (!z) return gc_NULL(av);
    1816      130403 :     z = FpX_to_mod(z, p);
    1817      130403 :     if (!pr || pr == ONLY_REM || pr == ONLY_DIVIDES)
    1818       71218 :       return gerepileupto(av, z);
    1819       59185 :     *pr = FpX_to_mod(*pr, p);
    1820       59185 :     return gc_all(av, 2, &z, pr);
    1821             :   }
    1822    19531910 :   switch(typ(y_lead))
    1823             :   {
    1824           0 :     case t_REAL:
    1825           0 :       y_lead = ginv(y_lead);
    1826           0 :       f = gmul; mod = NULL;
    1827           0 :       break;
    1828        2074 :     case t_INTMOD:
    1829        2074 :     case t_POLMOD: y_lead = ginv(y_lead);
    1830        2074 :       f = gmul; mod = gmodulo(gen_1, gel(y_lead,1));
    1831        2074 :       break;
    1832    19529836 :     default: if (gequal1(y_lead)) y_lead = NULL;
    1833    19529832 :       f = gdiv; mod = NULL;
    1834             :   }
    1835             : 
    1836    19531906 :   if (y_lead == NULL)
    1837    17518446 :     p2 = gel(x,dx+2);
    1838             :   else {
    1839             :     for(;;) {
    1840     2013460 :       p2 = f(gel(x,dx+2),y_lead);
    1841     2013459 :       p2 = simplify_shallow(p2);
    1842     2013459 :       if (!isexactzero(p2) || (--dx < 0)) break;
    1843             :     }
    1844     2013459 :     if (dx < dy) /* leading coeff of x was in fact zero */
    1845             :     {
    1846           0 :       if (pr == ONLY_DIVIDES) {
    1847           0 :         set_avma(av);
    1848           0 :         return (dx < 0)? pol_0(varn(x)) : NULL;
    1849             :       }
    1850           0 :       if (pr == ONLY_REM)
    1851             :       {
    1852           0 :         if (dx < 0)
    1853           0 :           return gerepilecopy(av, scalarpol(p2, varn(x)));
    1854             :         else
    1855             :         {
    1856             :           GEN t;
    1857           0 :           set_avma(av);
    1858           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1859           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1860           0 :           return t;
    1861             :         }
    1862             :       }
    1863           0 :       if (pr) /* cf ONLY_REM above */
    1864             :       {
    1865           0 :         if (dx < 0)
    1866             :         {
    1867           0 :           p2 = gclone(p2);
    1868           0 :           set_avma(av);
    1869           0 :           z = pol_0(varn(x));
    1870           0 :           x = scalarpol(p2, varn(x));
    1871           0 :           gunclone(p2);
    1872             :         }
    1873             :         else
    1874             :         {
    1875             :           GEN t;
    1876           0 :           set_avma(av);
    1877           0 :           z = pol_0(varn(x));
    1878           0 :           t = cgetg(dx + 3, t_POL); t[1] = x[1];
    1879           0 :           for (i = 2; i < dx + 3; i++) gel(t,i) = gcopy(gel(x,i));
    1880           0 :           x = t;
    1881             :         }
    1882           0 :         *pr = x;
    1883             :       }
    1884             :       else
    1885             :       {
    1886           0 :         set_avma(av);
    1887           0 :         z = pol_0(varn(x));
    1888             :       }
    1889           0 :       return z;
    1890             :     }
    1891             :   }
    1892             :   /* dx >= dy */
    1893    19531905 :   avy = avma;
    1894    19531905 :   dz = dx-dy;
    1895    19531905 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    1896    19531897 :   x += 2;
    1897    19531897 :   z += 2;
    1898    19531897 :   y += 2;
    1899    19531897 :   gel(z,dz) = gcopy(p2);
    1900             : 
    1901    54103054 :   for (i=dx-1; i>=dy; i--)
    1902             :   {
    1903    34571368 :     av1=avma; p1=gel(x,i);
    1904  1140512167 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1905    34566212 :     if (y_lead) p1 = simplify(f(p1,y_lead));
    1906             : 
    1907    34566212 :     if (isrationalzero(p1)) { set_avma(av1); p1 = gen_0; }
    1908             :     else
    1909    24575845 :       p1 = avma==av1? gcopy(p1): gerepileupto(av1,p1);
    1910    34570708 :     gel(z,i-dy) = p1;
    1911             :   }
    1912    19531686 :   if (!pr) return gerepileupto(av,z-2);
    1913             : 
    1914    12554583 :   rem = (GEN)avma; av1 = (pari_sp)new_chunk(dx+3);
    1915    12919370 :   for (sx=0; ; i--)
    1916             :   {
    1917    12919370 :     p1 = gel(x,i);
    1918             :     /* we always enter this loop at least once */
    1919    32443607 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1920    12918156 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1921    12918156 :     if (!gequal0(p1)) { sx = 1; break; } /* remainder is nonzero */
    1922     3588099 :     if (!isexactzero(p1)) break;
    1923     3437768 :     if (!i) break;
    1924      364730 :     set_avma(av1);
    1925             :   }
    1926    12553884 :   if (pr == ONLY_DIVIDES)
    1927             :   {
    1928        7931 :     if (sx) return gc_NULL(av);
    1929        7910 :     set_avma((pari_sp)rem); return gerepileupto(av,z-2);
    1930             :   }
    1931    12545953 :   lr=i+3; rem -= lr;
    1932    12545953 :   if (avma==av1) { set_avma((pari_sp)rem); p1 = gcopy(p1); }
    1933    12502348 :   else p1 = gerepileupto((pari_sp)rem,p1);
    1934    12546644 :   rem[0] = evaltyp(t_POL) | _evallg(lr);
    1935    12546644 :   rem[1] = z[-1];
    1936    12546644 :   rem += 2;
    1937    12546644 :   gel(rem,i) = p1;
    1938    17275535 :   for (i--; i>=0; i--)
    1939             :   {
    1940     4728897 :     av1=avma; p1 = gel(x,i);
    1941    13592013 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    1942     4728420 :     if (mod && avma==av1) p1 = gmul(p1,mod);
    1943     4728636 :     gel(rem,i) = avma==av1? gcopy(p1):gerepileupto(av1,p1);
    1944             :   }
    1945    12546638 :   rem -= 2;
    1946    12546638 :   if (!sx) (void)normalizepol_lg(rem, lr);
    1947    12546860 :   if (pr == ONLY_REM) return gerepileupto(av,rem);
    1948     6764206 :   z -= 2;
    1949             :   {
    1950     6764206 :     GEN *gptr[2]; gptr[0]=&z; gptr[1]=&rem;
    1951     6764206 :     gerepilemanysp(av,avy,gptr,2); *pr = rem; return z;
    1952             :   }
    1953             : }
    1954             : 
    1955             : GEN
    1956    14235034 : RgX_divrem(GEN x, GEN y, GEN *pr)
    1957             : {
    1958    14235034 :   if (pr == ONLY_REM) return RgX_rem(x, y);
    1959    14235034 :   return RgX_divrem_i(x, y, pr);
    1960             : }
    1961             : 
    1962             : /* x and y in (R[Y]/T)[X]  (lifted), T in R[Y]. y preferably monic */
    1963             : GEN
    1964      154815 : RgXQX_divrem(GEN x, GEN y, GEN T, GEN *pr)
    1965             : {
    1966             :   long vx, dx, dy, dz, i, j, sx, lr;
    1967             :   pari_sp av0, av, tetpil;
    1968             :   GEN z,p1,rem,lead;
    1969             : 
    1970      154815 :   if (!signe(y)) pari_err_INV("RgXQX_divrem",y);
    1971      154815 :   vx = varn(x);
    1972      154815 :   dx = degpol(x);
    1973      154815 :   dy = degpol(y);
    1974      154815 :   if (dx < dy)
    1975             :   {
    1976       40052 :     if (pr)
    1977             :     {
    1978       40024 :       av0 = avma; x = RgXQX_red(x, T);
    1979       40024 :       if (pr == ONLY_DIVIDES) { set_avma(av0); return signe(x)? NULL: gen_0; }
    1980       40017 :       if (pr == ONLY_REM) return x;
    1981           0 :       *pr = x;
    1982             :     }
    1983          28 :     return pol_0(vx);
    1984             :   }
    1985      114763 :   lead = leading_coeff(y);
    1986      114763 :   if (!dy) /* y is constant */
    1987             :   {
    1988         602 :     if (pr && pr != ONLY_DIVIDES)
    1989             :     {
    1990           0 :       if (pr == ONLY_REM) return pol_0(vx);
    1991           0 :       *pr = pol_0(vx);
    1992             :     }
    1993         602 :     if (gequal1(lead)) return RgX_copy(x);
    1994           0 :     av0 = avma; x = gmul(x, ginvmod(lead,T)); tetpil = avma;
    1995           0 :     return gerepile(av0,tetpil,RgXQX_red(x,T));
    1996             :   }
    1997      114161 :   av0 = avma; dz = dx-dy;
    1998      114161 :   lead = gequal1(lead)? NULL: gclone(ginvmod(lead,T));
    1999      114161 :   set_avma(av0);
    2000      114161 :   z = cgetg(dz+3,t_POL); z[1] = x[1];
    2001      114161 :   x += 2; y += 2; z += 2;
    2002             : 
    2003      114161 :   p1 = gel(x,dx); av = avma;
    2004      114161 :   gel(z,dz) = lead? gerepileupto(av, grem(gmul(p1,lead), T)): gcopy(p1);
    2005      580224 :   for (i=dx-1; i>=dy; i--)
    2006             :   {
    2007      466063 :     av=avma; p1=gel(x,i);
    2008     2414790 :     for (j=i-dy+1; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    2009      466029 :     if (lead) p1 = gmul(grem(p1, T), lead);
    2010      466027 :     tetpil=avma; gel(z,i-dy) = gerepile(av,tetpil, grem(p1, T));
    2011             :   }
    2012      114161 :   if (!pr) { guncloneNULL(lead); return z-2; }
    2013             : 
    2014      112887 :   rem = (GEN)avma; av = (pari_sp)new_chunk(dx+3);
    2015      189733 :   for (sx=0; ; i--)
    2016             :   {
    2017      189733 :     p1 = gel(x,i);
    2018      656520 :     for (j=0; j<=i && j<=dz; j++) p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    2019      189734 :     tetpil=avma; p1 = grem(p1, T); if (!gequal0(p1)) { sx = 1; break; }
    2020      104097 :     if (!i) break;
    2021       76846 :     set_avma(av);
    2022             :   }
    2023      112887 :   if (pr == ONLY_DIVIDES)
    2024             :   {
    2025       27246 :     guncloneNULL(lead);
    2026       27246 :     if (sx) return gc_NULL(av0);
    2027       24675 :     return gc_const((pari_sp)rem, z-2);
    2028             :   }
    2029       85641 :   lr=i+3; rem -= lr;
    2030       85641 :   rem[0] = evaltyp(t_POL) | _evallg(lr);
    2031       85641 :   rem[1] = z[-1];
    2032       85641 :   p1 = gerepile((pari_sp)rem,tetpil,p1);
    2033       85641 :   rem += 2; gel(rem,i) = p1;
    2034      187079 :   for (i--; i>=0; i--)
    2035             :   {
    2036      101438 :     av=avma; p1 = gel(x,i);
    2037      337694 :     for (j=0; j<=i && j<=dz; j++)
    2038      236256 :       p1 = gsub(p1, gmul(gel(z,j),gel(y,i-j)));
    2039      101438 :     tetpil=avma; gel(rem,i) = gerepile(av,tetpil, grem(p1, T));
    2040             :   }
    2041       85641 :   rem -= 2;
    2042       85641 :   guncloneNULL(lead);
    2043       85641 :   if (!sx) (void)normalizepol_lg(rem, lr);
    2044       85641 :   if (pr == ONLY_REM) return gerepileupto(av0,rem);
    2045         161 :   *pr = rem; return z-2;
    2046             : }
    2047             : 
    2048             : /*******************************************************************/
    2049             : /*                                                                 */
    2050             : /*                        PSEUDO-DIVISION                          */
    2051             : /*                                                                 */
    2052             : /*******************************************************************/
    2053             : INLINE GEN
    2054     4054864 : rem(GEN c, GEN T)
    2055             : {
    2056     4054864 :   if (T && typ(c) == t_POL && varn(c) == varn(T)) c = RgX_rem(c, T);
    2057     4054865 :   return c;
    2058             : }
    2059             : 
    2060             : /* x, y, are ZYX, lc(y) is an integer, T is a ZY */
    2061             : int
    2062       17331 : ZXQX_dvd(GEN x, GEN y, GEN T)
    2063             : {
    2064             :   long dx, dy, i, T_ismonic;
    2065       17331 :   pari_sp av = avma, av2;
    2066             :   GEN y_lead;
    2067             : 
    2068       17331 :   if (!signe(y)) pari_err_INV("ZXQX_dvd",y);
    2069       17331 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2070       17331 :   if (typ(y_lead) == t_POL) y_lead = gel(y_lead, 2); /* t_INT */
    2071             :   /* if monic, no point in using pseudo-division */
    2072       17331 :   if (gequal1(y_lead)) return signe(RgXQX_rem(x, y, T)) == 0;
    2073       14790 :   T_ismonic = gequal1(leading_coeff(T));
    2074       14790 :   dx = degpol(x);
    2075       14790 :   if (dx < dy) return !signe(x);
    2076       14790 :   (void)new_chunk(2);
    2077       14790 :   x = RgX_recip_i(x)+2;
    2078       14790 :   y = RgX_recip_i(y)+2;
    2079             :   /* pay attention to sparse divisors */
    2080       30081 :   for (i = 1; i <= dy; i++)
    2081       15291 :     if (!signe(gel(y,i))) gel(y,i) = NULL;
    2082       14790 :   av2 = avma;
    2083             :   for (;;)
    2084       73005 :   {
    2085       87795 :     GEN m, x0 = gel(x,0), y0 = y_lead, cx = content(x0);
    2086       87795 :     x0 = gneg(x0);
    2087       87795 :     m = gcdii(cx, y0);
    2088       87795 :     if (!equali1(m))
    2089             :     {
    2090       85293 :       x0 = gdiv(x0, m);
    2091       85293 :       y0 = diviiexact(y0, m);
    2092       85292 :       if (equali1(y0)) y0 = NULL;
    2093             :     }
    2094      179362 :     for (i=1; i<=dy; i++)
    2095             :     {
    2096       91567 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    2097       91567 :       if (gel(y,i)) c = gadd(c, gmul(x0,gel(y,i)));
    2098       91568 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    2099       91568 :       gel(x,i) = c;
    2100             :     }
    2101      688156 :     for (   ; i<=dx; i++)
    2102             :     {
    2103      600361 :       GEN c = gel(x,i); if (y0) c = gmul(y0, c);
    2104      600361 :       if (typ(c) == t_POL) c = T_ismonic ? ZX_rem(c, T): RgX_rem(c, T);
    2105      600361 :       gel(x,i) = c;
    2106             :     }
    2107      102341 :     do { x++; dx--; } while (dx >= 0 && !signe(gel(x,0)));
    2108       87795 :     if (dx < dy) break;
    2109       73005 :     if (gc_needed(av2,1))
    2110             :     {
    2111          28 :       if(DEBUGMEM>1) pari_warn(warnmem,"ZXQX_dvd dx = %ld >= %ld",dx,dy);
    2112          28 :       gerepilecoeffs(av2,x,dx+1);
    2113             :     }
    2114             :   }
    2115       14790 :   return gc_bool(av, dx < 0);
    2116             : }
    2117             : 
    2118             : /* T either NULL or a t_POL. */
    2119             : GEN
    2120      631105 : RgXQX_pseudorem(GEN x, GEN y, GEN T)
    2121             : {
    2122      631105 :   long vx = varn(x), dx, dy, dz, i, lx, p;
    2123      631105 :   pari_sp av = avma, av2;
    2124             :   GEN y_lead;
    2125             : 
    2126      631105 :   if (!signe(y)) pari_err_INV("RgXQX_pseudorem",y);
    2127      631105 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2128             :   /* if monic, no point in using pseudo-division */
    2129      631105 :   if (gequal1(y_lead)) return T? RgXQX_rem(x, y, T): RgX_rem(x, y);
    2130      581627 :   dx = degpol(x);
    2131      581627 :   if (dx < dy) return RgX_copy(x);
    2132      581620 :   (void)new_chunk(2);
    2133      581621 :   x = RgX_recip_i(x)+2;
    2134      581622 :   y = RgX_recip_i(y)+2;
    2135             :   /* pay attention to sparse divisors */
    2136     1959493 :   for (i = 1; i <= dy; i++)
    2137     1377871 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2138      581622 :   dz = dx-dy; p = dz+1;
    2139      581622 :   av2 = avma;
    2140             :   for (;;)
    2141             :   {
    2142     1250517 :     gel(x,0) = gneg(gel(x,0)); p--;
    2143     4204536 :     for (i=1; i<=dy; i++)
    2144             :     {
    2145     2954026 :       GEN c = gmul(y_lead, gel(x,i));
    2146     2953957 :       if (gel(y,i)) c = gadd(c, gmul(gel(x,0),gel(y,i)));
    2147     2954007 :       gel(x,i) = rem(c, T);
    2148             :     }
    2149     2189099 :     for (   ; i<=dx; i++)
    2150             :     {
    2151      938581 :       GEN c = gmul(y_lead, gel(x,i));
    2152      938572 :       gel(x,i) = rem(c, T);
    2153             :     }
    2154     1297047 :     do { x++; dx--; } while (dx >= 0 && gequal0(gel(x,0)));
    2155     1250516 :     if (dx < dy) break;
    2156      668895 :     if (gc_needed(av2,1))
    2157             :     {
    2158           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudorem dx = %ld >= %ld",dx,dy);
    2159           0 :       gerepilecoeffs(av2,x,dx+1);
    2160             :     }
    2161             :   }
    2162      581621 :   if (dx < 0) return pol_0(vx);
    2163      580438 :   lx = dx+3; x -= 2;
    2164      580438 :   x[0] = evaltyp(t_POL) | _evallg(lx);
    2165      580438 :   x[1] = evalsigne(1) | evalvarn(vx);
    2166      580438 :   x = RgX_recip_i(x);
    2167      580438 :   if (p)
    2168             :   { /* multiply by y[0]^p   [beware dummy vars from FpX_FpXY_resultant] */
    2169       28856 :     GEN t = y_lead;
    2170       28856 :     if (T && typ(t) == t_POL && varn(t) == varn(T))
    2171           0 :       t = RgXQ_powu(t, p, T);
    2172             :     else
    2173       28856 :       t = gpowgs(t, p);
    2174       89166 :     for (i=2; i<lx; i++)
    2175             :     {
    2176       60310 :       GEN c = gmul(gel(x,i), t);
    2177       60310 :       gel(x,i) = rem(c,T);
    2178             :     }
    2179       28856 :     if (!T) return gerepileupto(av, x);
    2180             :   }
    2181      551582 :   return gerepilecopy(av, x);
    2182             : }
    2183             : 
    2184             : GEN
    2185      631105 : RgX_pseudorem(GEN x, GEN y) { return RgXQX_pseudorem(x,y, NULL); }
    2186             : 
    2187             : /* Compute z,r s.t lc(y)^(dx-dy+1) x = z y + r */
    2188             : GEN
    2189       12098 : RgXQX_pseudodivrem(GEN x, GEN y, GEN T, GEN *ptr)
    2190             : {
    2191       12098 :   long vx = varn(x), dx, dy, dz, i, iz, lx, lz, p;
    2192       12098 :   pari_sp av = avma, av2;
    2193             :   GEN z, r, ypow, y_lead;
    2194             : 
    2195       12098 :   if (!signe(y)) pari_err_INV("RgXQX_pseudodivrem",y);
    2196       12098 :   dy = degpol(y); y_lead = gel(y,dy+2);
    2197       12098 :   if (gequal1(y_lead)) return T? RgXQX_divrem(x,y, T, ptr): RgX_divrem(x,y, ptr);
    2198        7603 :   dx = degpol(x);
    2199        7603 :   if (dx < dy) { *ptr = RgX_copy(x); return pol_0(vx); }
    2200        7603 :   if (dx == dy)
    2201             :   {
    2202          98 :     GEN x_lead = gel(x,lg(x)-1);
    2203          98 :     x = RgX_renormalize_lg(leafcopy(x), lg(x)-1);
    2204          98 :     y = RgX_renormalize_lg(leafcopy(y), lg(y)-1);
    2205          98 :     r = RgX_sub(RgX_Rg_mul(x, y_lead), RgX_Rg_mul(y, x_lead));
    2206          98 :     *ptr = gerepileupto(av, r); return scalarpol(x_lead, vx);
    2207             :   }
    2208        7505 :   (void)new_chunk(2);
    2209        7505 :   x = RgX_recip_i(x)+2;
    2210        7505 :   y = RgX_recip_i(y)+2;
    2211             :   /* pay attention to sparse divisors */
    2212       39084 :   for (i = 1; i <= dy; i++)
    2213       31579 :     if (isexactzero(gel(y,i))) gel(y,i) = NULL;
    2214        7505 :   dz = dx-dy; p = dz+1;
    2215        7505 :   lz = dz+3;
    2216        7505 :   z = cgetg(lz, t_POL);
    2217        7505 :   z[1] = evalsigne(1) | evalvarn(vx);
    2218       28808 :   for (i = 2; i < lz; i++) gel(z,i) = gen_0;
    2219        7505 :   ypow = new_chunk(dz+1);
    2220        7505 :   gel(ypow,0) = gen_1;
    2221        7505 :   gel(ypow,1) = y_lead;
    2222       13798 :   for (i=2; i<=dz; i++)
    2223             :   {
    2224        6293 :     GEN c = gmul(gel(ypow,i-1), y_lead);
    2225        6293 :     gel(ypow,i) = rem(c,T);
    2226             :   }
    2227        7505 :   av2 = avma;
    2228        7505 :   for (iz=2;;)
    2229             :   {
    2230       15731 :     p--;
    2231       15731 :     gel(z,iz++) = rem(gmul(gel(x,0), gel(ypow,p)), T);
    2232       70167 :     for (i=1; i<=dy; i++)
    2233             :     {
    2234       54436 :       GEN c = gmul(y_lead, gel(x,i));
    2235       54436 :       if (gel(y,i)) c = gsub(c, gmul(gel(x,0),gel(y,i)));
    2236       54436 :       gel(x,i) = rem(c, T);
    2237             :     }
    2238       41310 :     for (   ; i<=dx; i++)
    2239             :     {
    2240       25579 :       GEN c = gmul(y_lead, gel(x,i));
    2241       25579 :       gel(x,i) = rem(c,T);
    2242             :     }
    2243       15731 :     x++; dx--;
    2244       21303 :     while (dx >= dy && gequal0(gel(x,0))) { x++; dx--; iz++; }
    2245       15731 :     if (dx < dy) break;
    2246        8226 :     if (gc_needed(av2,1))
    2247             :     {
    2248           0 :       GEN X = x-2;
    2249           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgX_pseudodivrem dx=%ld >= %ld",dx,dy);
    2250           0 :       X[0] = evaltyp(t_POL)|_evallg(dx+3); X[1] = z[1]; /* hack */
    2251           0 :       gerepileall(av2,2, &X, &z); x = X+2;
    2252             :     }
    2253             :   }
    2254       14050 :   while (dx >= 0 && gequal0(gel(x,0))) { x++; dx--; }
    2255        7505 :   if (dx < 0)
    2256         182 :     x = pol_0(vx);
    2257             :   else
    2258             :   {
    2259        7323 :     lx = dx+3; x -= 2;
    2260        7323 :     x[0] = evaltyp(t_POL) | _evallg(lx);
    2261        7323 :     x[1] = evalsigne(1) | evalvarn(vx);
    2262        7323 :     x = RgX_recip_i(x);
    2263             :   }
    2264        7505 :   z = RgX_recip_i(z);
    2265        7505 :   r = x;
    2266        7505 :   if (p)
    2267             :   {
    2268        3339 :     GEN c = gel(ypow,p); r = RgX_Rg_mul(r, c);
    2269        3339 :     if (T && typ(c) == t_POL && varn(c) == varn(T)) r = RgXQX_red(r, T);
    2270             :   }
    2271        7505 :   *ptr = r; return gc_all(av, 2, &z, ptr);
    2272             : }
    2273             : GEN
    2274       11860 : RgX_pseudodivrem(GEN x, GEN y, GEN *ptr)
    2275       11860 : { return RgXQX_pseudodivrem(x,y,NULL,ptr); }
    2276             : 
    2277             : GEN
    2278           0 : RgXQX_mul(GEN x, GEN y, GEN T)
    2279           0 : { return RgXQX_red(RgX_mul(x,y), T); }
    2280             : GEN
    2281   733974052 : RgX_Rg_mul(GEN x, GEN y) { pari_APPLY_pol(gmul(y, gel(x,i))); }
    2282             : GEN
    2283      141351 : RgX_mul2n(GEN x, long n) { pari_APPLY_pol(gmul2n(gel(x,i), n)); }
    2284             : GEN
    2285       26425 : RgX_muls(GEN x, long y) { pari_APPLY_pol(gmulsg(y, gel(x,i))); }
    2286             : GEN
    2287          35 : RgXQX_RgXQ_mul(GEN x, GEN y, GEN T) { return RgXQX_red(RgX_Rg_mul(x,y), T); }
    2288             : GEN
    2289         133 : RgXQV_RgXQ_mul(GEN v, GEN x, GEN T) { return RgXQV_red(RgV_Rg_mul(v,x), T); }
    2290             : 
    2291             : GEN
    2292           0 : RgXQX_sqr(GEN x, GEN T) { return RgXQX_red(RgX_sqr(x), T); }
    2293             : 
    2294             : GEN
    2295           0 : RgXQX_powers(GEN P, long n, GEN T)
    2296             : {
    2297           0 :   GEN v = cgetg(n+2, t_VEC);
    2298             :   long i;
    2299           0 :   gel(v, 1) = pol_1(varn(T));
    2300           0 :   if (n==0) return v;
    2301           0 :   gel(v, 2) = gcopy(P);
    2302           0 :   for (i = 2; i <= n; i++) gel(v,i+1) = RgXQX_mul(P, gel(v,i), T);
    2303           0 :   return v;
    2304             : }
    2305             : 
    2306             : static GEN
    2307      623439 : _add(void *data, GEN x, GEN y) { (void)data; return RgX_add(x, y); }
    2308             : static GEN
    2309           0 : _sub(void *data, GEN x, GEN y) { (void)data; return RgX_sub(x, y); }
    2310             : static GEN
    2311       67310 : _sqr(void *data, GEN x) { return RgXQ_sqr(x, (GEN)data); }
    2312             : static GEN
    2313       88307 : _pow(void *data, GEN x, GEN n) { return RgXQ_pow(x, n, (GEN)data); }
    2314             : static GEN
    2315      276835 : _mul(void *data, GEN x, GEN y) { return RgXQ_mul(x,y, (GEN)data); }
    2316             : static GEN
    2317     1040579 : _cmul(void *data, GEN P, long a, GEN x) { (void)data; return RgX_Rg_mul(x,gel(P,a+2)); }
    2318             : static GEN
    2319     1019286 : _one(void *data) { return pol_1(varn((GEN)data)); }
    2320             : static GEN
    2321        1232 : _zero(void *data) { return pol_0(varn((GEN)data)); }
    2322             : static GEN
    2323      724481 : _red(void *data, GEN x) { (void)data; return gcopy(x); }
    2324             : 
    2325             : static struct bb_algebra RgXQ_algebra = { _red, _add, _sub,
    2326             :               _mul, _sqr, _one, _zero };
    2327             : 
    2328             : GEN
    2329           0 : RgX_RgXQV_eval(GEN Q, GEN x, GEN T)
    2330             : {
    2331           0 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)T,&RgXQ_algebra,_cmul);
    2332             : }
    2333             : 
    2334             : GEN
    2335      417397 : RgX_RgXQ_eval(GEN Q, GEN x, GEN T)
    2336             : {
    2337      417397 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2338      417398 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)T,&RgXQ_algebra,_cmul);
    2339             : }
    2340             : 
    2341             : /* mod X^n */
    2342             : struct modXn {
    2343             :   long v; /* varn(X) */
    2344             :   long n;
    2345             : } ;
    2346             : static GEN
    2347       11893 : _sqrXn(void *data, GEN x) {
    2348       11893 :   struct modXn *S = (struct modXn*)data;
    2349       11893 :   return RgXn_sqr(x, S->n);
    2350             : }
    2351             : static GEN
    2352        4528 : _mulXn(void *data, GEN x, GEN y) {
    2353        4528 :   struct modXn *S = (struct modXn*)data;
    2354        4528 :   return RgXn_mul(x,y, S->n);
    2355             : }
    2356             : static GEN
    2357        1939 : _oneXn(void *data) {
    2358        1939 :   struct modXn *S = (struct modXn*)data;
    2359        1939 :   return pol_1(S->v);
    2360             : }
    2361             : static GEN
    2362           0 : _zeroXn(void *data) {
    2363           0 :   struct modXn *S = (struct modXn*)data;
    2364           0 :   return pol_0(S->v);
    2365             : }
    2366             : static struct bb_algebra RgXn_algebra = { _red, _add, _sub, _mulXn, _sqrXn,
    2367             :                                           _oneXn, _zeroXn };
    2368             : 
    2369             : GEN
    2370         483 : RgXn_powers(GEN x, long m, long n)
    2371             : {
    2372         483 :   long d = degpol(x);
    2373         483 :   int use_sqr = (d<<1) >= n;
    2374             :   struct modXn S;
    2375         483 :   S.v = varn(x); S.n = n;
    2376         483 :   return gen_powers(x,m,use_sqr,(void*)&S,_sqrXn,_mulXn,_oneXn);
    2377             : }
    2378             : 
    2379             : GEN
    2380        2286 : RgXn_powu_i(GEN x, ulong m, long n)
    2381             : {
    2382             :   struct modXn S;
    2383             :   long v;
    2384        2286 :   if (n == 0) return x;
    2385        2286 :   v = RgX_valrem(x, &x);
    2386        2286 :   if (v) { n -= m * v; if (n <= 0) return pol_0(varn(x)); }
    2387        2265 :   S.v = varn(x); S.n = n;
    2388        2265 :   x = gen_powu_i(x, m, (void*)&S,_sqrXn,_mulXn);
    2389        2265 :   if (v) x = RgX_shift_shallow(x, m * v);
    2390        2265 :   return x;
    2391             : }
    2392             : GEN
    2393           0 : RgXn_powu(GEN x, ulong m, long n)
    2394             : {
    2395             :   pari_sp av;
    2396           0 :   if (n == 0) return gcopy(x);
    2397           0 :   av = avma; return gerepilecopy(av, RgXn_powu_i(x, m, n));
    2398             : }
    2399             : 
    2400             : GEN
    2401         966 : RgX_RgXnV_eval(GEN Q, GEN x, long n)
    2402             : {
    2403             :   struct modXn S;
    2404         966 :   S.v = varn(gel(x,2)); S.n = n;
    2405         966 :   return gen_bkeval_powers(Q,degpol(Q),x,(void*)&S,&RgXn_algebra,_cmul);
    2406             : }
    2407             : 
    2408             : GEN
    2409           0 : RgX_RgXn_eval(GEN Q, GEN x, long n)
    2410             : {
    2411           0 :   int use_sqr = 2*degpol(x) >= n;
    2412             :   struct modXn S;
    2413           0 :   S.v = varn(x); S.n = n;
    2414           0 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2415             : }
    2416             : 
    2417             : /* Q(x) mod t^n, x in R[t], n >= 1 */
    2418             : GEN
    2419        5159 : RgXn_eval(GEN Q, GEN x, long n)
    2420             : {
    2421        5159 :   long d = degpol(x);
    2422             :   int use_sqr;
    2423             :   struct modXn S;
    2424        5159 :   if (d == 1 && isrationalzero(gel(x,2)))
    2425             :   {
    2426        5152 :     GEN y = RgX_unscale(Q, gel(x,3));
    2427        5152 :     setvarn(y, varn(x)); return y;
    2428             :   }
    2429           7 :   S.v = varn(x);
    2430           7 :   S.n = n;
    2431           7 :   use_sqr = (d<<1) >= n;
    2432           7 :   return gen_bkeval(Q,degpol(Q),x,use_sqr,(void*)&S,&RgXn_algebra,_cmul);
    2433             : }
    2434             : 
    2435             : /* (f*g mod t^n) \ t^n2, assuming 2*n2 >= n */
    2436             : static GEN
    2437     2614505 : RgXn_mulhigh(GEN f, GEN g, long n2, long n)
    2438             : {
    2439     2614505 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2440     2614505 :   return RgX_add(RgX_mulhigh_i(fl, g, n2), RgXn_mul(fh, g, n - n2));
    2441             : }
    2442             : 
    2443             : /* (f^2 mod t^n) \ t^n2, assuming 2*n2 >= n */
    2444             : static GEN
    2445          14 : RgXn_sqrhigh(GEN f, long n2, long n)
    2446             : {
    2447          14 :   GEN F = RgX_blocks(f, n2, 2), fl = gel(F,1), fh = gel(F,2);
    2448          14 :   return RgX_add(RgX_mulhigh_i(fl, f, n2), RgXn_mul(fh, f, n - n2));
    2449             : }
    2450             : 
    2451             : static GEN
    2452     2188922 : RgXn_div_gen(GEN g, GEN f, long e)
    2453             : {
    2454             :   pari_sp av;
    2455             :   ulong mask;
    2456             :   GEN W, a;
    2457     2188922 :   long v = varn(f), n = 1;
    2458             : 
    2459     2188922 :   if (!signe(f)) pari_err_INV("RgXn_inv",f);
    2460     2188915 :   a = ginv(gel(f,2));
    2461     2188912 :   if (e == 1 && !g) return scalarpol(a, v);
    2462     2184026 :   else if (e == 2 && !g)
    2463             :   {
    2464             :     GEN b;
    2465      836930 :     if (degpol(f) <= 0 || gequal0(b = gel(f,3))) return scalarpol(a, v);
    2466      246250 :     b = gneg(b);
    2467      246252 :     if (!gequal1(a)) b = gmul(b, gsqr(a));
    2468      246251 :     return deg1pol(b, a, v);
    2469             :   }
    2470     1347096 :   av = avma;
    2471     1347096 :   W = scalarpol_shallow(a,v);
    2472     1347096 :   mask = quadratic_prec_mask(e);
    2473     3955623 :   while (mask > 1)
    2474             :   {
    2475             :     GEN u, fr;
    2476     2608527 :     long n2 = n;
    2477     2608527 :     n<<=1; if (mask & 1) n--;
    2478     2608527 :     mask >>= 1;
    2479     2608527 :     fr = RgXn_red_shallow(f, n);
    2480     2608527 :     if (mask>1 || !g)
    2481             :     {
    2482     1324547 :       u = RgXn_mul(W, RgXn_mulhigh(fr, W, n2, n), n-n2);
    2483     1324547 :       W = RgX_sub(W, RgX_shift_shallow(u, n2));
    2484             :     }
    2485             :     else
    2486             :     {
    2487     1283980 :       GEN y = RgXn_mul(g, W, n), yt =  RgXn_red_shallow(y, n-n2);
    2488     1283980 :       u = RgXn_mul(yt, RgXn_mulhigh(fr,  W, n2, n), n-n2);
    2489     1283980 :       W = RgX_sub(y, RgX_shift_shallow(u, n2));
    2490             :     }
    2491     2608527 :     if (gc_needed(av,2))
    2492             :     {
    2493           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_inv, e = %ld", n);
    2494           0 :       W = gerepileupto(av, W);
    2495             :     }
    2496             :   }
    2497     1347096 :   return W;
    2498             : }
    2499             : 
    2500             : static GEN
    2501         119 : RgXn_div_FpX(GEN x, GEN y, long e, GEN p)
    2502             : {
    2503             :   GEN r;
    2504         119 :   if (lgefint(p) == 3)
    2505             :   {
    2506         119 :     ulong pp = uel(p, 2);
    2507         119 :     if (pp == 2)
    2508          14 :       r = F2x_to_ZX(F2xn_div(RgX_to_F2x(x), RgX_to_F2x(y), e));
    2509             :     else
    2510         105 :       r = Flx_to_ZX_inplace(Flxn_div(RgX_to_Flx(x, pp), RgX_to_Flx(y, pp), e, pp));
    2511             :   }
    2512             :   else
    2513           0 :     r = FpXn_div(RgX_to_FpX(x, p), RgX_to_FpX(y, p), e, p);
    2514         119 :   return FpX_to_mod(r, p);
    2515             : }
    2516             : 
    2517             : static GEN
    2518           7 : RgXn_div_FpXQX(GEN x, GEN y, long n, GEN pol, GEN p)
    2519             : {
    2520           7 :   GEN r, T = RgX_to_FpX(pol, p);
    2521           7 :   if (signe(T) == 0) pari_err_OP("/", x, y);
    2522           7 :   r = FpXQXn_div(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), n, T, p);
    2523           7 :   return FpXQX_to_mod(r, T, p);
    2524             : }
    2525             : 
    2526             : static GEN
    2527          91 : RgXn_inv_FpX(GEN x, long e, GEN p)
    2528             : {
    2529             :   GEN r;
    2530          91 :   if (lgefint(p) == 3)
    2531             :   {
    2532          91 :     ulong pp = uel(p, 2);
    2533          91 :     if (pp == 2)
    2534          28 :       r = F2x_to_ZX(F2xn_inv(RgX_to_F2x(x), e));
    2535             :     else
    2536          63 :       r = Flx_to_ZX_inplace(Flxn_inv(RgX_to_Flx(x, pp), e, pp));
    2537             :   }
    2538             :   else
    2539           0 :     r = FpXn_inv(RgX_to_FpX(x, p), e, p);
    2540          91 :   return FpX_to_mod(r, p);
    2541             : }
    2542             : 
    2543             : static GEN
    2544           0 : RgXn_inv_FpXQX(GEN x, long n, GEN pol, GEN p)
    2545             : {
    2546           0 :   GEN r, T = RgX_to_FpX(pol, p);
    2547           0 :   if (signe(T) == 0) pari_err_OP("/", gen_1, x);
    2548           0 :   r = FpXQXn_inv(RgX_to_FpXQX(x, T, p), n, T, p);
    2549           0 :   return FpXQX_to_mod(r, T, p);
    2550             : }
    2551             : 
    2552             : static GEN
    2553      905033 : RgXn_inv_fast(GEN x, long e)
    2554             : {
    2555             :   GEN p, pol;
    2556             :   long pa;
    2557      905033 :   long t = RgX_type(x,&p,&pol,&pa);
    2558      905034 :   switch(t)
    2559             :   {
    2560          91 :     case t_INTMOD: return RgXn_inv_FpX(x, e, p);
    2561           0 :     case RgX_type_code(t_POLMOD, t_INTMOD):
    2562           0 :                    return RgXn_inv_FpXQX(x, e, pol, p);
    2563      904943 :     default:       return NULL;
    2564             :   }
    2565             : }
    2566             : 
    2567             : static GEN
    2568     1284106 : RgXn_div_fast(GEN x, GEN y, long e)
    2569             : {
    2570             :   GEN p, pol;
    2571             :   long pa;
    2572     1284106 :   long t = RgX_type2(x,y,&p,&pol,&pa);
    2573     1284106 :   switch(t)
    2574             :   {
    2575         119 :     case t_INTMOD: return RgXn_div_FpX(x, y, e, p);
    2576           7 :     case RgX_type_code(t_POLMOD, t_INTMOD):
    2577           7 :                    return RgXn_div_FpXQX(x, y, e, pol, p);
    2578     1283980 :     default:       return NULL;
    2579             :   }
    2580             : }
    2581             : 
    2582             : GEN
    2583     1284106 : RgXn_div_i(GEN g, GEN f, long e)
    2584             : {
    2585     1284106 :   GEN h = RgXn_div_fast(g, f, e);
    2586     1284106 :   if (h) return h;
    2587     1283980 :   return RgXn_div_gen(g, f, e);
    2588             : }
    2589             : 
    2590             : GEN
    2591      584443 : RgXn_div(GEN g, GEN f, long e)
    2592             : {
    2593      584443 :   pari_sp av = avma;
    2594      584443 :   return gerepileupto(av, RgXn_div_i(g, f, e));
    2595             : }
    2596             : 
    2597             : GEN
    2598      905033 : RgXn_inv_i(GEN f, long e)
    2599             : {
    2600      905033 :   GEN h = RgXn_inv_fast(f, e);
    2601      905034 :   if (h) return h;
    2602      904943 :   return RgXn_div_gen(NULL, f, e);
    2603             : }
    2604             : 
    2605             : GEN
    2606      807450 : RgXn_inv(GEN f, long e)
    2607             : {
    2608      807450 :   pari_sp av = avma;
    2609      807450 :   return gerepileupto(av, RgXn_inv_i(f, e));
    2610             : }
    2611             : 
    2612             : /* intformal(x^n*S) / x^(n+1) */
    2613             : static GEN
    2614       56417 : RgX_integXn(GEN x, long n)
    2615      114435 : { pari_APPLY_pol_normalized(gdivgs(gel(x,i), n+i-1)); }
    2616             : 
    2617             : GEN
    2618       52941 : RgXn_expint(GEN h, long e)
    2619             : {
    2620       52941 :   pari_sp av = avma, av2;
    2621       52941 :   long v = varn(h), n;
    2622       52941 :   GEN f = pol_1(v), g;
    2623             :   ulong mask;
    2624             : 
    2625       52939 :   if (!signe(h)) return f;
    2626       50448 :   g = pol_1(v);
    2627       50447 :   n = 1; mask = quadratic_prec_mask(e);
    2628       50447 :   av2 = avma;
    2629       56418 :   for (;mask>1;)
    2630             :   {
    2631             :     GEN u, w;
    2632       56418 :     long n2 = n;
    2633       56418 :     n<<=1; if (mask & 1) n--;
    2634       56418 :     mask >>= 1;
    2635       56418 :     u = RgXn_mul(g, RgX_mulhigh_i(f, RgXn_red_shallow(h, n2-1), n2-1), n-n2);
    2636       56417 :     u = RgX_add(u, RgX_shift_shallow(RgXn_red_shallow(h, n-1), 1-n2));
    2637       56417 :     w = RgXn_mul(f, RgX_integXn(u, n2-1), n-n2);
    2638       56416 :     f = RgX_add(f, RgX_shift_shallow(w, n2));
    2639       56416 :     if (mask<=1) break;
    2640        5971 :     u = RgXn_mul(g, RgXn_mulhigh(f, g, n2, n), n-n2);
    2641        5971 :     g = RgX_sub(g, RgX_shift_shallow(u, n2));
    2642        5971 :     if (gc_needed(av2,2))
    2643             :     {
    2644           0 :       if (DEBUGMEM>1) pari_warn(warnmem,"RgXn_expint, e = %ld", n);
    2645           0 :       gerepileall(av2, 2, &f, &g);
    2646             :     }
    2647             :   }
    2648       50445 :   return gerepileupto(av, f);
    2649             : }
    2650             : 
    2651             : GEN
    2652           0 : RgXn_exp(GEN h, long e)
    2653             : {
    2654           0 :   long d = degpol(h);
    2655           0 :   if (d < 0) return pol_1(varn(h));
    2656           0 :   if (!d || !gequal0(gel(h,2)))
    2657           0 :     pari_err_DOMAIN("RgXn_exp","valuation", "<", gen_1, h);
    2658           0 :   return RgXn_expint(RgX_deriv(h), e);
    2659             : }
    2660             : 
    2661             : GEN
    2662         154 : RgXn_reverse(GEN f, long e)
    2663             : {
    2664         154 :   pari_sp av = avma, av2;
    2665             :   ulong mask;
    2666             :   GEN fi, a, df, W, an;
    2667         154 :   long v = varn(f), n=1;
    2668         154 :   if (degpol(f)<1 || !gequal0(gel(f,2)))
    2669           0 :     pari_err_INV("serreverse",f);
    2670         154 :   fi = ginv(gel(f,3));
    2671         154 :   a = deg1pol_shallow(fi,gen_0,v);
    2672         154 :   if (e <= 2) return gerepilecopy(av, a);
    2673         133 :   W = scalarpol(fi,v);
    2674         133 :   df = RgX_deriv(f);
    2675         133 :   mask = quadratic_prec_mask(e);
    2676         133 :   av2 = avma;
    2677         616 :   for (;mask>1;)
    2678             :   {
    2679             :     GEN u, fa, fr;
    2680         483 :     long n2 = n, rt;
    2681         483 :     n<<=1; if (mask & 1) n--;
    2682         483 :     mask >>= 1;
    2683         483 :     fr = RgXn_red_shallow(f, n);
    2684         483 :     rt = brent_kung_optpow(degpol(fr), 4, 3);
    2685         483 :     an = RgXn_powers(a, rt, n);
    2686         483 :     if (n>1)
    2687             :     {
    2688         483 :       long n4 = (n2+1)>>1;
    2689         483 :       GEN dfr = RgXn_red_shallow(df, n2);
    2690         483 :       dfr = RgX_RgXnV_eval(dfr, RgXnV_red_shallow(an, n2), n2);
    2691         483 :       u = RgX_shift(RgX_Rg_sub(RgXn_mul(W, dfr, n2), gen_1), -n4);
    2692         483 :       W = RgX_sub(W, RgX_shift(RgXn_mul(u, W, n2-n4), n4));
    2693             :     }
    2694         483 :     fa = RgX_sub(RgX_RgXnV_eval(fr, an, n), pol_x(v));
    2695         483 :     fa = RgX_shift(fa, -n2);
    2696         483 :     a = RgX_sub(a, RgX_shift(RgXn_mul(W, fa, n-n2), n2));
    2697         483 :     if (gc_needed(av2,2))
    2698             :     {
    2699           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_reverse, e = %ld", n);
    2700           0 :       gerepileall(av2, 2, &a, &W);
    2701             :     }
    2702             :   }
    2703         133 :   return gerepileupto(av, a);
    2704             : }
    2705             : 
    2706             : GEN
    2707           7 : RgXn_sqrt(GEN h, long e)
    2708             : {
    2709           7 :   pari_sp av = avma, av2;
    2710           7 :   long v = varn(h), n = 1;
    2711           7 :   GEN f = scalarpol(gen_1, v), df = f;
    2712           7 :   ulong mask = quadratic_prec_mask(e);
    2713           7 :   if (degpol(h)<0 || !gequal1(gel(h,2)))
    2714           0 :     pari_err_SQRTN("RgXn_sqrt",h);
    2715           7 :   av2 = avma;
    2716             :   while(1)
    2717           7 :   {
    2718          14 :     long n2 = n, m;
    2719             :     GEN g;
    2720          14 :     n<<=1; if (mask & 1) n--;
    2721          14 :     mask >>= 1;
    2722          14 :     m = n-n2;
    2723          14 :     g = RgX_sub(RgXn_sqrhigh(f, n2, n), RgX_shift_shallow(RgXn_red_shallow(h, n),-n2));
    2724          14 :     f = RgX_sub(f, RgX_shift_shallow(RgXn_mul(gmul2n(df, -1), g, m), n2));
    2725          14 :     if (mask==1) return gerepileupto(av, f);
    2726           7 :     g = RgXn_mul(df, RgXn_mulhigh(df, f, n2, n), m);
    2727           7 :     df = RgX_sub(df, RgX_shift_shallow(g, n2));
    2728           7 :     if (gc_needed(av2,2))
    2729             :     {
    2730           0 :       if(DEBUGMEM>1) pari_warn(warnmem,"RgXn_sqrt, e = %ld", n);
    2731           0 :       gerepileall(av2, 2, &f, &df);
    2732             :     }
    2733             :   }
    2734             : }
    2735             : 
    2736             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2737             : GEN
    2738      116836 : RgXQ_powu(GEN x, ulong n, GEN T)
    2739             : {
    2740      116836 :   pari_sp av = avma;
    2741             : 
    2742      116836 :   if (!n) return pol_1(varn(x));
    2743       68533 :   if (n == 1) return RgX_copy(x);
    2744       19845 :   x = gen_powu_i(x, n, (void*)T, &_sqr, &_mul);
    2745       19845 :   return gerepilecopy(av, x);
    2746             : }
    2747             : /* x,T in Rg[X], n in N, compute lift(x^n mod T)) */
    2748             : GEN
    2749      102160 : RgXQ_pow(GEN x, GEN n, GEN T)
    2750             : {
    2751             :   pari_sp av;
    2752      102160 :   long s = signe(n);
    2753             : 
    2754      102160 :   if (!s) return pol_1(varn(x));
    2755      102160 :   if (is_pm1(n) == 1)
    2756       88307 :     return (s < 0)? RgXQ_inv(x, T): RgX_copy(x);
    2757       13853 :   av = avma;
    2758       13853 :   if (s < 0) x = RgXQ_inv(x, T);
    2759       13853 :   x = gen_pow_i(x, n, (void*)T, &_sqr, &_mul);
    2760       13853 :   return gerepilecopy(av, x);
    2761             : }
    2762             : static GEN
    2763      200590 : _ZXQsqr(void *data, GEN x) { return ZXQ_sqr(x, (GEN)data); }
    2764             : static GEN
    2765      111189 : _ZXQmul(void *data, GEN x, GEN y) { return ZXQ_mul(x,y, (GEN)data); }
    2766             : 
    2767             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2768             : GEN
    2769       11781 : ZXQ_powers(GEN x, long l, GEN T)
    2770             : {
    2771       11781 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2772       11781 :   return gen_powers(x, l, use_sqr, (void *)T,_ZXQsqr,_ZXQmul,_one);
    2773             : }
    2774             : 
    2775             : /* x,T in Z[X], n in N, compute lift(x^n mod T)) */
    2776             : GEN
    2777      168504 : ZXQ_powu(GEN x, ulong n, GEN T)
    2778             : {
    2779      168504 :   pari_sp av = avma;
    2780             : 
    2781      168504 :   if (!n) return pol_1(varn(x));
    2782      168504 :   if (n == 1) return ZX_copy(x);
    2783      113868 :   x = gen_powu_i(x, n, (void*)T, &_ZXQsqr, &_ZXQmul);
    2784      113874 :   return gerepilecopy(av, x);
    2785             : }
    2786             : 
    2787             : /* generates the list of powers of x of degree 0,1,2,...,l*/
    2788             : GEN
    2789        8575 : RgXQ_powers(GEN x, long l, GEN T)
    2790             : {
    2791        8575 :   int use_sqr = 2*degpol(x) >= degpol(T);
    2792        8575 :   return gen_powers(x, l, use_sqr, (void *)T,_sqr,_mul,_one);
    2793             : }
    2794             : 
    2795             : GEN
    2796        7118 : RgXQV_factorback(GEN L, GEN e, GEN T)
    2797             : {
    2798        7118 :   return gen_factorback(L, e, (void*)T, &_mul, &_pow, &_one);
    2799             : }
    2800             : 
    2801             : /* a in K = Q[X]/(T), returns [a^0, ..., a^n] */
    2802             : GEN
    2803        9226 : QXQ_powers(GEN a, long n, GEN T)
    2804             : {
    2805             :   GEN den, v;
    2806        9226 :   if (!isint1(leading_coeff(T))) return RgXQ_powers(a, n, T);
    2807        9198 :   v = ZXQ_powers(Q_remove_denom(a, &den), n, T);
    2808             :   /* den*a integral; v[i+1] = (den*a)^i in K */
    2809        9198 :   if (den)
    2810             :   { /* restore denominators */
    2811        5605 :     GEN d = den;
    2812             :     long i;
    2813        5605 :     gel(v,2) = a;
    2814       27237 :     for (i=3; i<=n+1; i++) {
    2815       21632 :       d = mulii(d,den);
    2816       21632 :       gel(v,i) = RgX_Rg_div(gel(v,i), d);
    2817             :     }
    2818             :   }
    2819        9198 :   return v;
    2820             : }
    2821             : 
    2822             : static GEN
    2823        3234 : do_QXQ_eval(GEN v, long imin, GEN a, GEN T)
    2824             : {
    2825        3234 :   long l, i, m = 0;
    2826             :   GEN dz, z;
    2827        3234 :   GEN V = cgetg_copy(v, &l);
    2828       10843 :   for (i = imin; i < l; i++)
    2829             :   {
    2830        7609 :     GEN c = gel(v, i);
    2831        7609 :     if (typ(c) == t_POL) m = maxss(m, degpol(c));
    2832             :   }
    2833        3234 :   z = Q_remove_denom(QXQ_powers(a, m, T), &dz);
    2834        3465 :   for (i = 1; i < imin; i++) V[i] = v[i];
    2835       10843 :   for (i = imin; i < l; i++)
    2836             :   {
    2837        7609 :     GEN c = gel(v,i);
    2838        7609 :     if (typ(c) == t_POL) c = QX_ZXQV_eval(c, z, dz);
    2839        7609 :     gel(V,i) = c;
    2840             :   }
    2841        3234 :   return V;
    2842             : }
    2843             : /* [ s(a mod T) | s <- lift(v) ], a,T are QX, v a QXV */
    2844             : GEN
    2845        3003 : QXV_QXQ_eval(GEN v, GEN a, GEN T)
    2846        3003 : { return do_QXQ_eval(v, 1, a, T); }
    2847             : 
    2848             : GEN
    2849         231 : QXY_QXQ_evalx(GEN v, GEN a, GEN T)
    2850         231 : { return normalizepol(do_QXQ_eval(v, 2, a, T)); }
    2851             : 
    2852             : GEN
    2853        2961 : RgXQ_matrix_pow(GEN y, long n, long m, GEN P)
    2854             : {
    2855        2961 :   return RgXV_to_RgM(RgXQ_powers(y,m-1,P),n);
    2856             : }
    2857             : 
    2858             : GEN
    2859        5144 : RgXQ_norm(GEN x, GEN T)
    2860             : {
    2861             :   pari_sp av;
    2862        5144 :   long dx = degpol(x);
    2863             :   GEN L, y;
    2864        5144 :   if (degpol(T)==0) return gpowgs(x,0);
    2865        5137 :   av = avma; y = resultant(T, x);
    2866        5137 :   L = leading_coeff(T);
    2867        5137 :   if (gequal1(L) || !signe(x)) return y;
    2868           0 :   return gerepileupto(av, gdiv(y, gpowgs(L, dx)));
    2869             : }
    2870             : 
    2871             : GEN
    2872         476 : RgXQ_trace(GEN x, GEN T)
    2873             : {
    2874         476 :   pari_sp av = avma;
    2875             :   GEN dT, z;
    2876             :   long n;
    2877         476 :   if (degpol(T)==0) return gmulgs(x,0);
    2878         469 :   dT = RgX_deriv(T); n = degpol(dT);
    2879         469 :   z = RgXQ_mul(x, dT, T);
    2880         469 :   if (degpol(z)<n) return gc_const(av, gen_0);
    2881         420 :   return gerepileupto(av, gdiv(gel(z,2+n), gel(T,3+n)));
    2882             : }
    2883             : 
    2884             : GEN
    2885     2788190 : RgX_blocks(GEN P, long n, long m)
    2886             : {
    2887     2788190 :   GEN z = cgetg(m+1,t_VEC);
    2888     2788190 :   long i,j, k=2, l = lg(P);
    2889     8558148 :   for(i=1; i<=m; i++)
    2890             :   {
    2891     5769958 :     GEN zi = cgetg(n+2,t_POL);
    2892     5769958 :     zi[1] = P[1];
    2893     5769958 :     gel(z,i) = zi;
    2894    16997926 :     for(j=2; j<n+2; j++)
    2895    11227968 :       gel(zi, j) = k==l ? gen_0 : gel(P,k++);
    2896     5769958 :     zi = RgX_renormalize_lg(zi, n+2);
    2897             :   }
    2898     2788190 :   return z;
    2899             : }
    2900             : 
    2901             : /* write p(X) = e(X^2) + Xo(X^2), shallow function */
    2902             : void
    2903    49759161 : RgX_even_odd(GEN p, GEN *pe, GEN *po)
    2904             : {
    2905    49759161 :   long n = degpol(p), v = varn(p), n0, n1, i;
    2906             :   GEN p0, p1;
    2907             : 
    2908    49759601 :   if (n <= 0) { *pe = RgX_copy(p); *po = zeropol(v); return; }
    2909             : 
    2910    49666324 :   n0 = (n>>1)+1; n1 = n+1 - n0; /* n1 <= n0 <= n1+1 */
    2911    49666324 :   p0 = cgetg(n0+2, t_POL); p0[1] = evalvarn(v)|evalsigne(1);
    2912    49668620 :   p1 = cgetg(n1+2, t_POL); p1[1] = evalvarn(v)|evalsigne(1);
    2913   146794653 :   for (i=0; i<n1; i++)
    2914             :   {
    2915    97124242 :     p0[2+i] = p[2+(i<<1)];
    2916    97124242 :     p1[2+i] = p[3+(i<<1)];
    2917             :   }
    2918    49670411 :   if (n1 != n0)
    2919    18008105 :     p0[2+i] = p[2+(i<<1)];
    2920    49670411 :   *pe = normalizepol(p0);
    2921    49674927 :   *po = normalizepol(p1);
    2922             : }
    2923             : 
    2924             : /* write p(X) = a_0(X^k) + Xa_1(X^k) + ... + X^(k-1)a_{k-1}(X^k), shallow function */
    2925             : GEN
    2926       43631 : RgX_splitting(GEN p, long k)
    2927             : {
    2928       43631 :   long n = degpol(p), v = varn(p), m, i, j, l;
    2929             :   GEN r;
    2930             : 
    2931       43631 :   m = n/k;
    2932       43631 :   r = cgetg(k+1,t_VEC);
    2933      234059 :   for(i=1; i<=k; i++)
    2934             :   {
    2935      190428 :     gel(r,i) = cgetg(m+3, t_POL);
    2936      190428 :     mael(r,i,1) = evalvarn(v)|evalsigne(1);
    2937             :   }
    2938      571991 :   for (j=1, i=0, l=2; i<=n; i++)
    2939             :   {
    2940      528360 :     gmael(r,j,l) = gel(p,2+i);
    2941      528360 :     if (j==k) { j=1; l++; } else j++;
    2942             :   }
    2943      234059 :   for(i=1; i<=k; i++)
    2944      190428 :     gel(r,i) = normalizepol_lg(gel(r,i),i<j?l+1:l);
    2945       43631 :   return r;
    2946             : }
    2947             : 
    2948             : /*******************************************************************/
    2949             : /*                                                                 */
    2950             : /*                        Kronecker form                           */
    2951             : /*                                                                 */
    2952             : /*******************************************************************/
    2953             : 
    2954             : /* z in R[Y] representing an elt in R[X,Y] mod T(Y) in Kronecker form,
    2955             :  * i.e subst(lift(z), x, y^(2deg(z)-1)). Recover the "real" z, with
    2956             :  * normalized coefficients */
    2957             : GEN
    2958      187337 : Kronecker_to_mod(GEN z, GEN T)
    2959             : {
    2960      187337 :   long i,j,lx,l = lg(z), N = (degpol(T)<<1) + 1;
    2961      187337 :   GEN x, t = cgetg(N,t_POL);
    2962      187337 :   t[1] = T[1];
    2963      187337 :   lx = (l-2) / (N-2); x = cgetg(lx+3,t_POL);
    2964      187337 :   x[1] = z[1];
    2965      187337 :   T = RgX_copy(T);
    2966     1406637 :   for (i=2; i<lx+2; i++, z+= N-2)
    2967             :   {
    2968     5793872 :     for (j=2; j<N; j++) gel(t,j) = gel(z,j);
    2969     1219300 :     gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2970             :   }
    2971      187337 :   N = (l-2) % (N-2) + 2;
    2972      475782 :   for (j=2; j<N; j++) t[j] = z[j];
    2973      187337 :   gel(x,i) = mkpolmod(RgX_rem(normalizepol_lg(t,N), T), T);
    2974      187337 :   return normalizepol_lg(x, i+1);
    2975             : }
    2976             : 
    2977             : /*******************************************************************/
    2978             : /*                                                                 */
    2979             : /*                        Domain detection                         */
    2980             : /*                                                                 */
    2981             : /*******************************************************************/
    2982             : 
    2983             : static GEN
    2984      666063 : RgX_mul_FpX(GEN x, GEN y, GEN p)
    2985             : {
    2986      666063 :   pari_sp av = avma;
    2987             :   GEN r;
    2988      666063 :   if (lgefint(p) == 3)
    2989             :   {
    2990      614707 :     ulong pp = uel(p, 2);
    2991      614707 :     r = Flx_to_ZX_inplace(Flx_mul(RgX_to_Flx(x, pp),
    2992             :                                   RgX_to_Flx(y, pp), pp));
    2993             :   }
    2994             :   else
    2995       51356 :     r = FpX_mul(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    2996      666063 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    2997      561230 :   return gerepileupto(av, FpX_to_mod(r, p));
    2998             : }
    2999             : 
    3000             : static GEN
    3001         616 : zero_FpXQX_mod(GEN pol, GEN p, long v)
    3002             : {
    3003         616 :   GEN r = cgetg(3,t_POL);
    3004         616 :   r[1] = evalvarn(v);
    3005         616 :   gel(r,2) = mkpolmod(mkintmod(gen_0, icopy(p)), gcopy(pol));
    3006         616 :   return r;
    3007             : }
    3008             : 
    3009             : static GEN
    3010        2079 : RgX_mul_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    3011             : {
    3012        2079 :   pari_sp av = avma;
    3013             :   long dT;
    3014             :   GEN kx, ky, r;
    3015        2079 :   GEN T = RgX_to_FpX(pol, p);
    3016        2079 :   if (signe(T)==0) pari_err_OP("*", x, y);
    3017        2072 :   dT = degpol(T);
    3018        2072 :   kx = RgXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    3019        2072 :   ky = RgXX_to_Kronecker(RgX_to_FpXQX(y, T, p), dT);
    3020        2072 :   r = FpX_mul(kx, ky, p);
    3021        2072 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    3022        1463 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    3023             : }
    3024             : 
    3025             : static GEN
    3026      455653 : RgX_liftred(GEN x, GEN T)
    3027      455653 : { return RgXQX_red(liftpol_shallow(x), T); }
    3028             : 
    3029             : static GEN
    3030      172260 : RgX_mul_QXQX(GEN x, GEN y, GEN T)
    3031             : {
    3032      172260 :   pari_sp av = avma;
    3033      172260 :   long dT = degpol(T);
    3034      172260 :   GEN r = QX_mul(RgXX_to_Kronecker(RgX_liftred(x, T), dT),
    3035             :                  RgXX_to_Kronecker(RgX_liftred(y, T), dT));
    3036      172260 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    3037             : }
    3038             : 
    3039             : static GEN
    3040        1421 : RgX_sqr_FpX(GEN x, GEN p)
    3041             : {
    3042        1421 :   pari_sp av = avma;
    3043             :   GEN r;
    3044        1421 :   if (lgefint(p) == 3)
    3045             :   {
    3046        1217 :     ulong pp = uel(p, 2);
    3047        1217 :     r = Flx_to_ZX_inplace(Flx_sqr(RgX_to_Flx(x, pp), pp));
    3048             :   }
    3049             :   else
    3050         204 :     r = FpX_sqr(RgX_to_FpX(x, p), p);
    3051        1421 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    3052        1281 :   return gerepileupto(av, FpX_to_mod(r, p));
    3053             : }
    3054             : 
    3055             : static GEN
    3056         196 : RgX_sqr_FpXQX(GEN x, GEN pol, GEN p)
    3057             : {
    3058         196 :   pari_sp av = avma;
    3059             :   long dT;
    3060         196 :   GEN kx, r, T = RgX_to_FpX(pol, p);
    3061         196 :   if (signe(T)==0) pari_err_OP("*",x,x);
    3062         189 :   dT = degpol(T);
    3063         189 :   kx = RgXX_to_Kronecker(RgX_to_FpXQX(x, T, p), dT);
    3064         189 :   r = FpX_sqr(kx, p);
    3065         189 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    3066         189 :   return gerepileupto(av, Kronecker_to_mod(FpX_to_mod(r, p), pol));
    3067             : }
    3068             : 
    3069             : static GEN
    3070       13425 : RgX_sqr_QXQX(GEN x, GEN T)
    3071             : {
    3072       13425 :   pari_sp av = avma;
    3073       13425 :   long dT = degpol(T);
    3074       13425 :   GEN r = QX_sqr(RgXX_to_Kronecker(RgX_liftred(x, T), dT));
    3075       13425 :   return gerepileupto(av, Kronecker_to_mod(r, T));
    3076             : }
    3077             : 
    3078             : static GEN
    3079       68320 : RgX_rem_FpX(GEN x, GEN y, GEN p)
    3080             : {
    3081       68320 :   pari_sp av = avma;
    3082             :   GEN r;
    3083       68320 :   if (lgefint(p) == 3)
    3084             :   {
    3085       51684 :     ulong pp = uel(p, 2);
    3086       51684 :     r = Flx_to_ZX_inplace(Flx_rem(RgX_to_Flx(x, pp),
    3087             :                                   RgX_to_Flx(y, pp), pp));
    3088             :   }
    3089             :   else
    3090       16636 :     r = FpX_rem(RgX_to_FpX(x, p), RgX_to_FpX(y, p), p);
    3091       68320 :   if (signe(r)==0) { set_avma(av); return zero_FpX_mod(p, varn(x)); }
    3092       31024 :   return gerepileupto(av, FpX_to_mod(r, p));
    3093             : }
    3094             : 
    3095             : static GEN
    3096       48854 : RgX_rem_QXQX(GEN x, GEN y, GEN T)
    3097             : {
    3098       48854 :   pari_sp av = avma;
    3099             :   GEN r;
    3100       48854 :   r = RgXQX_rem(RgX_liftred(x, T), RgX_liftred(y, T), T);
    3101       48854 :   return gerepilecopy(av, QXQX_to_mod_shallow(r, T));
    3102             : }
    3103             : static GEN
    3104          70 : RgX_rem_FpXQX(GEN x, GEN y, GEN pol, GEN p)
    3105             : {
    3106          70 :   pari_sp av = avma;
    3107             :   GEN r;
    3108          70 :   GEN T = RgX_to_FpX(pol, p);
    3109          70 :   if (signe(T) == 0) pari_err_OP("%", x, y);
    3110          63 :   if (lgefint(p) == 3)
    3111             :   {
    3112          55 :     ulong pp = uel(p, 2);
    3113          55 :     GEN Tp = ZX_to_Flx(T, pp);
    3114          55 :     r = FlxX_to_ZXX(FlxqX_rem(RgX_to_FlxqX(x, Tp, pp),
    3115             :                               RgX_to_FlxqX(y, Tp, pp), Tp, pp));
    3116             :   }
    3117             :   else
    3118           8 :     r = FpXQX_rem(RgX_to_FpXQX(x, T, p), RgX_to_FpXQX(y, T, p), T, p);
    3119          63 :   if (signe(r)==0) { set_avma(av); return zero_FpXQX_mod(pol, p, varn(x)); }
    3120          56 :   return gerepileupto(av, FpXQX_to_mod(r, T, p));
    3121             : }
    3122             : 
    3123             : static GEN
    3124   116839748 : RgX_mul_fast(GEN x, GEN y)
    3125             : {
    3126             :   GEN p, pol;
    3127             :   long pa;
    3128   116839748 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    3129   116840013 :   switch(t)
    3130             :   {
    3131    64859533 :     case t_INT:    return ZX_mul(x,y);
    3132     3355255 :     case t_FRAC:   return QX_mul(x,y);
    3133      105176 :     case t_FFELT:  return FFX_mul(x, y, pol);
    3134      666063 :     case t_INTMOD: return RgX_mul_FpX(x, y, p);
    3135      172260 :     case RgX_type_code(t_POLMOD, t_INT):
    3136             :     case RgX_type_code(t_POLMOD, t_FRAC):
    3137      172260 :                    return RgX_mul_QXQX(x, y, pol);
    3138        2038 :     case RgX_type_code(t_POLMOD, t_INTMOD):
    3139        2038 :                    return RgX_mul_FpXQX(x, y, pol, p);
    3140    47679688 :     default:       return NULL;
    3141             :   }
    3142             : }
    3143             : static GEN
    3144     1386759 : RgX_sqr_fast(GEN x)
    3145             : {
    3146             :   GEN p, pol;
    3147             :   long pa;
    3148     1386759 :   long t = RgX_type(x,&p,&pol,&pa);
    3149     1386764 :   switch(t)
    3150             :   {
    3151     1250056 :     case t_INT:    return ZX_sqr(x);
    3152       81581 :     case t_FRAC:   return QX_sqr(x);
    3153        2499 :     case t_FFELT:  return FFX_sqr(x, pol);
    3154        1421 :     case t_INTMOD: return RgX_sqr_FpX(x, p);
    3155       13425 :     case RgX_type_code(t_POLMOD, t_INT):
    3156             :     case RgX_type_code(t_POLMOD, t_FRAC):
    3157       13425 :                    return RgX_sqr_QXQX(x, pol);
    3158         192 :     case RgX_type_code(t_POLMOD, t_INTMOD):
    3159         192 :                    return RgX_sqr_FpXQX(x, pol, p);
    3160       37590 :     default:       return NULL;
    3161             :   }
    3162             : }
    3163             : 
    3164             : static GEN
    3165    14647026 : RgX_rem_fast(GEN x, GEN y)
    3166             : {
    3167             :   GEN p, pol;
    3168             :   long pa;
    3169    14647026 :   long t = RgX_type2(x,y, &p,&pol,&pa);
    3170    14647139 :   switch(t)
    3171             :   {
    3172     3424460 :     case t_INT:    return ZX_is_monic(y) ? ZX_rem(x,y): NULL;
    3173     1837620 :     case t_FRAC:   return RgX_is_ZX(y) && ZX_is_monic(y) ? QX_ZX_rem(x,y): NULL;
    3174          84 :     case t_FFELT:  return FFX_rem(x, y, pol);
    3175       68320 :     case t_INTMOD: return RgX_rem_FpX(x, y, p);
    3176       48854 :     case RgX_type_code(t_POLMOD, t_INT):
    3177             :     case RgX_type_code(t_POLMOD, t_FRAC):
    3178       48854 :                    return RgX_rem_QXQX(x, y, pol);
    3179          64 :     case RgX_type_code(t_POLMOD, t_INTMOD):
    3180          64 :                    return RgX_rem_FpXQX(x, y, pol, p);
    3181     9267737 :     default:       return NULL;
    3182             :   }
    3183             : }
    3184             : 
    3185             : GEN
    3186   103107446 : RgX_mul(GEN x, GEN y)
    3187             : {
    3188   103107446 :   GEN z = RgX_mul_fast(x,y);
    3189   103107575 :   if (!z) z = RgX_mul_i(x,y);
    3190   103107146 :   return z;
    3191             : }
    3192             : 
    3193             : GEN
    3194     1374607 : RgX_sqr(GEN x)
    3195             : {
    3196     1374607 :   GEN z = RgX_sqr_fast(x);
    3197     1374605 :   if (!z) z = RgX_sqr_i(x);
    3198     1374605 :   return z;
    3199             : }
    3200             : 
    3201             : GEN
    3202    14647031 : RgX_rem(GEN x, GEN y)
    3203             : {
    3204    14647031 :   GEN z = RgX_rem_fast(x, y);
    3205    14647125 :   if (!z) z = RgX_divrem_i(x, y, ONLY_REM);
    3206    14647085 :   return z;
    3207             : }
    3208             : 
    3209             : static GEN
    3210           0 : _RgX_mul(void* E, GEN x, GEN y)
    3211           0 : { (void) E; return RgX_mul(x, y); }
    3212             : 
    3213             : GEN
    3214           0 : RgXV_prod(GEN V)
    3215           0 : { return gen_product(V, NULL, &_RgX_mul); }
    3216             : 
    3217             : GEN
    3218    11061361 : RgXn_mul(GEN f, GEN g, long n)
    3219             : {
    3220    11061361 :   pari_sp av = avma;
    3221    11061361 :   GEN h = RgX_mul_fast(f,g);
    3222    11061364 :   if (!h) return RgXn_mul2(f,g,n);
    3223     1783548 :   if (degpol(h) < n) return h;
    3224      382941 :   return gerepilecopy(av, RgXn_red_shallow(h, n));
    3225             : }
    3226             : 
    3227             : GEN
    3228       12152 : RgXn_sqr(GEN f, long n)
    3229             : {
    3230       12152 :   pari_sp av = avma;
    3231       12152 :   GEN g = RgX_sqr_fast(f);
    3232       12152 :   if (!g) return RgXn_sqr2(f,n);
    3233       11837 :   if (degpol(g) < n) return g;
    3234       10640 :   return gerepilecopy(av, RgXn_red_shallow(g, n));
    3235             : }
    3236             : 
    3237             : /* (f*g) \/ x^n */
    3238             : GEN
    3239     2670937 : RgX_mulhigh_i(GEN f, GEN g, long n)
    3240             : {
    3241     2670937 :   GEN h = RgX_mul_fast(f,g);
    3242     2670937 :   return h? RgX_shift_shallow(h, -n): RgX_mulhigh_i2(f,g,n);
    3243             : }
    3244             : 
    3245             : /* (f*g) \/ x^n */
    3246             : GEN
    3247           0 : RgX_sqrhigh_i(GEN f, long n)
    3248             : {
    3249           0 :   GEN h = RgX_sqr_fast(f);
    3250           0 :   return h? RgX_shift_shallow(h, -n): RgX_sqrhigh_i2(f,n);
    3251             : }

Generated by: LCOV version 1.16