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 - gen1.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.18.1 lcov report (development 30005-fc14bb602a) Lines: 1884 2008 93.8 %
Date: 2025-02-18 09:22:46 Functions: 103 103 100.0 %
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             : /********************************************************************/
      16             : /**                                                                **/
      17             : /**                      GENERIC OPERATIONS                        **/
      18             : /**                         (first part)                           **/
      19             : /**                                                                **/
      20             : /********************************************************************/
      21             : #include "pari.h"
      22             : #include "paripriv.h"
      23             : 
      24             : #define DEBUGLEVEL DEBUGLEVEL_mod
      25             : 
      26             : /* assume z[1] was created last */
      27             : #define fix_frac_if_int(z) if (equali1(gel(z,2)))\
      28             :   z = gerepileupto((pari_sp)(z+3), gel(z,1));
      29             : 
      30             : /* assume z[1] was created last */
      31             : #define fix_frac_if_int_GC(z,tetpil) { if (equali1(gel(z,2)))\
      32             :   z = gerepileupto((pari_sp)(z+3), gel(z,1));\
      33             : else\
      34             :   gerepilecoeffssp((pari_sp)z, tetpil, z+1, 2); }
      35             : 
      36             : static void
      37         105 : warn_coercion(GEN x, GEN y, GEN z)
      38             : {
      39         105 :   if (DEBUGLEVEL)
      40          56 :    pari_warn(warner,"coercing quotient rings; moduli %Ps and %Ps -> %Ps",x,y,z);
      41         105 : }
      42             : 
      43             : static long
      44          28 : kro_quad(GEN x, GEN y)
      45          28 : { pari_sp av=avma; return gc_long(av, kronecker(quad_disc(x), y)); }
      46             : 
      47             : /* is -1 not a square in Zp, assume p prime */
      48             : INLINE int
      49          42 : Zp_nosquare_m1(GEN p) { return (mod4(p) & 2); /* 2 or 3 mod 4 */ }
      50             : 
      51             : static GEN addsub_pp(GEN x, GEN y, GEN(*op)(GEN,GEN));
      52             : static GEN mulpp(GEN x, GEN y);
      53             : static GEN divpp(GEN x, GEN y);
      54             : /* Argument codes for inline routines
      55             :  * c: complex, p: padic, q: quad, f: floating point (REAL, some complex)
      56             :  * R: without imaginary part (INT, REAL, INTMOD, FRAC, PADIC if -1 not square)
      57             :  * T: some type (to be converted to PADIC)
      58             :  */
      59             : static GEN
      60   304016307 : addRc(GEN x, GEN y) {
      61   304016307 :   GEN z = cgetg(3,t_COMPLEX);
      62   303996017 :   gel(z,1) = gadd(x,gel(y,1));
      63   303996782 :   gel(z,2) = gcopy(gel(y,2)); return z;
      64             : }
      65             : static GEN
      66   351578685 : mulRc(GEN x, GEN y) {
      67   351578685 :   GEN z = cgetg(3,t_COMPLEX);
      68   351531037 :   gel(z,1) = isintzero(gel(y,1))? gen_0: gmul(x,gel(y,1));
      69   351938867 :   gel(z,2) = gmul(x,gel(y,2)); return z;
      70             : }
      71             : /* for INTMODs: can't simplify when Re(x) = gen_0 */
      72             : static GEN
      73          49 : mulRc_direct(GEN x, GEN y) {
      74          49 :   GEN z = cgetg(3,t_COMPLEX);
      75          49 :   gel(z,1) = gmul(x,gel(y,1));
      76          49 :   gel(z,2) = gmul(x,gel(y,2)); return z;
      77             : }
      78             : static GEN
      79     2359618 : divRc(GEN x, GEN y) {
      80     2359618 :   GEN t = gdiv(x, cxnorm(y)), mt = gneg(t); /* left on stack for efficiency */
      81     2359615 :   GEN z = cgetg(3,t_COMPLEX);
      82     2359619 :   gel(z,1) = isintzero(gel(y,1))? gen_0: gmul(t, gel(y,1));
      83     2359608 :   gel(z,2) = gmul(mt, gel(y,2));
      84     2359608 :   return z;
      85             : }
      86             : static GEN
      87    23127364 : divcR(GEN x, GEN y) {
      88    23127364 :   GEN z = cgetg(3,t_COMPLEX);
      89    23127284 :   gel(z,1) = isintzero(gel(x,1))? gen_0: gdiv(gel(x,1), y);
      90    23126878 :   gel(z,2) = gdiv(gel(x,2), y); return z;
      91             : }
      92             : static GEN
      93        1274 : addRq(GEN x, GEN y) {
      94        1274 :   GEN z = cgetg(4,t_QUAD);
      95        1274 :   gel(z,1) = ZX_copy(gel(y,1));
      96        1274 :   gel(z,2) = gadd(x, gel(y,2));
      97        1274 :   gel(z,3) = gcopy(gel(y,3)); return z;
      98             : }
      99             : static GEN
     100        2310 : mulRq(GEN x, GEN y) {
     101        2310 :   GEN z = cgetg(4,t_QUAD);
     102        2310 :   gel(z,1) = ZX_copy(gel(y,1));
     103        2310 :   gel(z,2) = gmul(x,gel(y,2));
     104        2310 :   gel(z,3) = gmul(x,gel(y,3)); return z;
     105             : }
     106             : static GEN
     107          77 : addqf(GEN x, GEN y, long prec) { pari_sp av = avma;
     108          77 :   long i = gexpo(x) - gexpo(y);
     109          77 :   if (i > 0) prec += nbits2extraprec(i);
     110          77 :   return gerepileupto(av, gadd(y, quadtofp(x, prec)));
     111             : }
     112             : static GEN
     113    37513010 : mulrfrac(GEN x, GEN y)
     114             : {
     115             :   pari_sp av;
     116    37513010 :   GEN z, a = gel(y,1), b = gel(y,2);
     117    37513010 :   if (is_pm1(a)) /* frequent special case */
     118             :   {
     119    12657331 :     z = divri(x, b); if (signe(a) < 0) togglesign(z);
     120    12657884 :     return z;
     121             :   }
     122    24855754 :   av = avma;
     123    24855754 :   return gerepileuptoleaf(av, divri(mulri(x,a), b));
     124             : }
     125             : static GEN
     126          28 : mulqf(GEN x, GEN y, long prec) { pari_sp av = avma;
     127          28 :   return gerepileupto(av, gmul(y, quadtofp(x, prec)));
     128             : }
     129             : static GEN
     130          63 : divqf(GEN x, GEN y, long prec) { pari_sp av = avma;
     131          63 :   return gerepileupto(av, gdiv(quadtofp(x,prec), y));
     132             : }
     133             : static GEN
     134          42 : divfq(GEN x, GEN y, long prec) { pari_sp av = avma;
     135          42 :   return gerepileupto(av, gdiv(x, quadtofp(y,prec)));
     136             : }
     137             : /* y PADIC, x + y by converting x to padic */
     138             : static GEN
     139           7 : addTp(GEN x, GEN y) { pari_sp av = avma; GEN z;
     140             : 
     141           7 :   if (!valp(y)) z = cvtop2(x,y);
     142             :   else {
     143           7 :     long l = signe(padic_u(y))? valp(y) + precp(y): valp(y);
     144           7 :     z  = cvtop(x, padic_p(y), l);
     145             :   }
     146           7 :   return gerepileupto(av, addsub_pp(z, y, addii));
     147             : }
     148             : /* y PADIC, x * y by converting x to padic */
     149             : static GEN
     150     4982898 : mulTp(GEN x, GEN y) { pari_sp av = avma;
     151     4982898 :   return gerepileupto(av, mulpp(cvtop2(x,y), y));
     152             : }
     153             : /* y PADIC, non zero x / y by converting x to padic */
     154             : static GEN
     155        3660 : divTp(GEN x, GEN y) { pari_sp av = avma;
     156        3660 :   return gerepileupto(av, divpp(cvtop2(x,y), y));
     157             : }
     158             : /* x PADIC, x / y by converting y to padic. Assume x != 0; otherwise y
     159             :  * converted to O(p^e) and division by 0 */
     160             : static GEN
     161     1204910 : divpT(GEN x, GEN y) { pari_sp av = avma;
     162     1204910 :   return gerepileupto(av, divpp(x, cvtop2(y,x)));
     163             : }
     164             : 
     165             : /* z := Mod(x,X) + Mod(y,X) [ t_INTMOD preallocated ], x,y,X INT, 0 <= x,y < X
     166             :  * clean memory from z on */
     167             : static GEN
     168     3550175 : add_intmod_same(GEN z, GEN X, GEN x, GEN y) {
     169     3550175 :   if (lgefint(X) == 3) {
     170     3227045 :     ulong u = Fl_add(itou(x),itou(y), X[2]);
     171     3227045 :     set_avma((pari_sp)z); gel(z,2) = utoi(u);
     172             :   }
     173             :   else {
     174      323130 :     GEN u = addii(x,y); if (cmpii(u, X) >= 0) u = subii(u, X);
     175      323130 :     gel(z,2) = gerepileuptoint((pari_sp)z, u);
     176             :   }
     177     3550175 :   gel(z,1) = icopy(X); return z;
     178             : }
     179             : static GEN
     180     1053721 : sub_intmod_same(GEN z, GEN X, GEN x, GEN y) {
     181     1053721 :   if (lgefint(X) == 3) {
     182      679687 :     ulong u = Fl_sub(itou(x),itou(y), X[2]);
     183      679687 :     set_avma((pari_sp)z); gel(z,2) = utoi(u);
     184             :   }
     185             :   else {
     186      374034 :     GEN u = subii(x,y); if (signe(u) < 0) u = addii(u, X);
     187      374034 :     gel(z,2) = gerepileuptoint((pari_sp)z, u);
     188             :   }
     189     1053721 :   gel(z,1) = icopy(X); return z;
     190             : }
     191             : /* cf add_intmod_same */
     192             : static GEN
     193     2835705 : mul_intmod_same(GEN z, GEN X, GEN x, GEN y) {
     194     2835705 :   if (lgefint(X) == 3) {
     195     2694320 :     ulong u = Fl_mul(itou(x),itou(y), X[2]);
     196     2694320 :     set_avma((pari_sp)z); gel(z,2) = utoi(u);
     197             :   }
     198             :   else
     199      141385 :     gel(z,2) = gerepileuptoint((pari_sp)z, remii(mulii(x,y), X) );
     200     2835704 :   gel(z,1) = icopy(X); return z;
     201             : }
     202             : /* cf add_intmod_same */
     203             : static GEN
     204      341964 : div_intmod_same(GEN z, GEN X, GEN x, GEN y)
     205             : {
     206      341964 :   if (lgefint(X) == 3) {
     207      319417 :     ulong m = uel(X,2), u = Fl_div(itou(x), itou(y), m);
     208      319354 :     set_avma((pari_sp)z); gel(z,2) = utoi(u);
     209             :   }
     210             :   else
     211       22547 :     gel(z,2) = gerepileuptoint((pari_sp)z, remii(mulii(x, Fp_inv(y,X)), X) );
     212      341900 :   gel(z,1) = icopy(X); return z;
     213             : }
     214             : 
     215             : /*******************************************************************/
     216             : /*                                                                 */
     217             : /*        REDUCTION to IRREDUCIBLE TERMS (t_FRAC/t_RFRAC)          */
     218             : /*                                                                 */
     219             : /* (static routines are not memory clean, but OK for gerepileupto) */
     220             : /*******************************************************************/
     221             : /* Compute the denominator of (1/y) * (n/d) = n/yd, y a "scalar".
     222             :  * Sanity check : avoid (1/2) / (Mod(1,2)*x + 1) "=" 1 / (0 * x + 1) */
     223             : static GEN
     224     9857748 : rfrac_denom_mul_scal(GEN d, GEN y)
     225             : {
     226     9857748 :   GEN D = RgX_Rg_mul(d, y);
     227     9857748 :   if (lg(D) != lg(d))
     228             :   { /* try to generate a meaningful diagnostic */
     229           0 :     D = gdiv(leading_coeff(d), y); /* should fail */
     230           0 :     pari_err_INV("gred_rfrac", y); /* better than nothing */
     231             :   }
     232     9857748 :   return D;
     233             : }
     234             : 
     235             : static int
     236    57698003 : Leading_is_neg(GEN x)
     237             : {
     238   122016343 :   while (typ(x) == t_POL) x = leading_coeff(x);
     239    57698003 :   return is_real_t(typ(x))? gsigne(x) < 0: 0;
     240             : }
     241             : 
     242             : static int
     243   154190368 : transtype(GEN x) { return x != gen_1 && typ(x) != t_PADIC; }
     244             : 
     245             : /* d a t_POL, n a coprime t_POL of same var or "scalar". Not memory clean */
     246             : GEN
     247    57702239 : gred_rfrac_simple(GEN n, GEN d)
     248             : {
     249             :   GEN _1n, _1d, c, cn, cd, z;
     250    57702239 :   long dd = degpol(d);
     251             : 
     252    57702239 :   if (dd <= 0)
     253             :   {
     254        4236 :     if (dd < 0) pari_err_INV("gred_rfrac_simple", d);
     255        4236 :     n = gdiv(n, gel(d,2));
     256        4236 :     if (typ(n) != t_POL || varn(n) != varn(d)) n = scalarpol(n, varn(d));
     257        4236 :     return n;
     258             :   }
     259    57698003 :   if (Leading_is_neg(d)) { d = gneg(d); n = gneg(n); }
     260    57698003 :   _1n = Rg_get_1(n);
     261    57698003 :   _1d = Rg_get_1(d);
     262    57698003 :   if (transtype(_1n) && !gidentical(_1n, _1d)) d = gmul(d, _1n);
     263    57698003 :   if (transtype(_1d) && !gidentical(_1n, _1d)) n = gmul(n, _1d);
     264    57697996 :   cd = content(d);
     265    59564600 :   while (typ(n) == t_POL && !degpol(n)) n = gel(n,2);
     266    57697996 :   cn = (typ(n) == t_POL && varn(n) == varn(d))? content(n): n;
     267    57697996 :   if (!gequal1(cd)) {
     268     6567981 :     d = RgX_Rg_div(d,cd);
     269     6567981 :     if (!gequal1(cn))
     270             :     {
     271     1294189 :       if (gequal0(cn)) {
     272          49 :         if (isexactzero(cn)) return scalarpol(cn, varn(d));
     273           0 :         n = (cn != n)? RgX_Rg_div(n,cd): gdiv(n, cd);
     274           0 :         c = gen_1;
     275             :       } else {
     276     1294140 :         n = (cn != n)? RgX_Rg_div(n,cn): gen_1;
     277     1294140 :         c = gdiv(cn,cd);
     278             :       }
     279             :     }
     280             :     else
     281     5273792 :       c = ginv(cd);
     282             :   } else {
     283    51130015 :     if (!gequal1(cn))
     284             :     {
     285     3277620 :       if (gequal0(cn)) {
     286         987 :         if (isexactzero(cn)) return scalarpol(cn, varn(d));
     287          21 :         c = gen_1;
     288             :       } else {
     289     3276633 :         n = (cn != n)? RgX_Rg_div(n,cn): gen_1;
     290     3276633 :         c = cn;
     291             :       }
     292             :     } else {
     293    47852395 :       GEN y = cgetg(3,t_RFRAC);
     294    47852395 :       gel(y,1) = gcopy(n);
     295    47852395 :       gel(y,2) = RgX_copy(d); return y;
     296             :     }
     297             :   }
     298             : 
     299     9844586 :   if (typ(c) == t_POL)
     300             :   {
     301      901227 :     z = c;
     302      939657 :     do { z = content(z); } while (typ(z) == t_POL);
     303      901227 :     cd = denom_i(z);
     304      901227 :     cn = gmul(c, cd);
     305             :   }
     306             :   else
     307             :   {
     308     8943359 :     cn = numer_i(c);
     309     8943359 :     cd = denom_i(c);
     310             :   }
     311     9844586 :   z = cgetg(3,t_RFRAC);
     312     9844586 :   gel(z,1) = gmul(n, cn);
     313     9844586 :   gel(z,2) = d = rfrac_denom_mul_scal(d, cd);
     314             :   /* can happen: Pol(O(17^-1)) / Pol([Mod(9,23), O(23^-3)]) */
     315     9844586 :   if (!signe(d)) pari_err_INV("gred_rfrac_simple", d);
     316     9844586 :   return z;
     317             : }
     318             : 
     319             : /* in rare cases x may be a t_POL, after 0/x for instance -> pol_0() */
     320             : static GEN
     321      207774 : fix_rfrac(GEN x, long d)
     322             : {
     323             :   GEN z, N, D;
     324      207774 :   if (!d || typ(x) == t_POL) return x;
     325      165872 :   z = cgetg(3, t_RFRAC);
     326      165872 :   N = gel(x,1);
     327      165872 :   D = gel(x,2);
     328      165872 :   if (d > 0) {
     329        6265 :     gel(z, 1) = (typ(N)==t_POL && varn(N)==varn(D))? RgX_shift(N,d)
     330      160447 :                                                    : monomialcopy(N,d,varn(D));
     331      160384 :     gel(z, 2) = RgX_copy(D);
     332             :   } else {
     333        5488 :     gel(z, 1) = gcopy(N);
     334        5488 :     gel(z, 2) = RgX_shift(D, -d);
     335             :   }
     336      165872 :   return z;
     337             : }
     338             : 
     339             : /* assume d != 0 */
     340             : static GEN
     341    44760431 : gred_rfrac2(GEN n, GEN d)
     342             : {
     343             :   GEN y, z, _1n, _1d;
     344             :   long v, vd, vn;
     345             : 
     346    44760431 :   n = simplify_shallow(n);
     347    44760431 :   if (isintzero(n)) return scalarpol(Rg_get_0(d), varn(d));
     348    37919527 :   d = simplify_shallow(d);
     349    37919527 :   if (typ(d) != t_POL) return gdiv(n,d);
     350    36740783 :   vd = varn(d);
     351    36740783 :   if (typ(n) != t_POL)
     352             :   {
     353    20620712 :     if (varncmp(vd, gvar(n)) >= 0) return gdiv(n,d);
     354    20619305 :     if (varncmp(vd, gvar2(n)) < 0) return gred_rfrac_simple(n,d);
     355           0 :     pari_err_BUG("gred_rfrac2 [incompatible variables]");
     356             :   }
     357    16120071 :   vn = varn(n);
     358    16120071 :   if (varncmp(vd, vn) < 0) return gred_rfrac_simple(n,d);
     359    15982362 :   if (varncmp(vd, vn) > 0) return RgX_Rg_div(n,d);
     360    15822589 :   _1n = Rg_get_1(n);
     361    15822589 :   _1d = Rg_get_1(d);
     362    15822589 :   if (transtype(_1n) && !gidentical(_1n, _1d)) d = gmul(d, _1n);
     363    15822589 :   if (transtype(_1d) && !gidentical(_1n, _1d)) n = gmul(n, _1d);
     364             : 
     365             :   /* now n and d are t_POLs in the same variable */
     366    15822589 :   v = RgX_valrem(n, &n) - RgX_valrem(d, &d);
     367    15822589 :   if (!degpol(d))
     368             :   {
     369    12784887 :     n = RgX_Rg_div(n,gel(d,2));
     370    12784887 :     return v? RgX_mulXn(n,v): n;
     371             :   }
     372             : 
     373             :   /* X does not divide gcd(n,d), deg(d) > 0 */
     374     3037702 :   if (!isinexact(n) && !isinexact(d))
     375             :   {
     376     3037478 :     y = RgX_divrem(n, d, &z);
     377     3037478 :     if (!signe(z)) { cgiv(z); return v? RgX_mulXn(y, v): y; }
     378      207550 :     z = RgX_gcd(d, z);
     379      207550 :     if (degpol(z)) { n = RgX_div(n,z); d = RgX_div(d,z); }
     380             :   }
     381      207774 :   return fix_rfrac(gred_rfrac_simple(n,d), v);
     382             : }
     383             : 
     384             : /* x,y t_INT, return x/y in reduced form */
     385             : GEN
     386   132157158 : Qdivii(GEN x, GEN y)
     387             : {
     388   132157158 :   pari_sp av = avma;
     389             :   GEN r, q;
     390             : 
     391   132157158 :   if (lgefint(y) == 3)
     392             :   {
     393   114172922 :     q = Qdiviu(x, y[2]);
     394   114171081 :     if (signe(y) > 0) return q;
     395    11013628 :     if (typ(q) == t_INT) togglesign(q); else togglesign_safe(&gel(q,1));
     396    11013862 :     return q;
     397             :   }
     398    17984236 :   if (is_pm1(y)) return (signe(y) < 0)? negi(x): icopy(x);
     399    17985733 :   if (equali1(x))
     400             :   {
     401     5182062 :     if (!signe(y)) pari_err_INV("gdiv",y);
     402     5181957 :     retmkfrac(signe(y) < 0? gen_m1: gen_1, absi(y));
     403             :   }
     404    12803711 :   q = dvmdii(x,y,&r);
     405    12803626 :   if (r == gen_0) return q; /* gen_0 intended */
     406     8378112 :   r = gcdii(y, r);
     407     8378124 :   if (lgefint(r) == 3)
     408             :   {
     409     7494518 :     ulong t = r[2];
     410     7494518 :     set_avma(av);
     411     7494485 :     if (t == 1) q = mkfraccopy(x,y);
     412             :     else
     413             :     {
     414     2650118 :       q = cgetg(3,t_FRAC);
     415     2650254 :       gel(q,1) = diviuexact(x,t);
     416     2650216 :       gel(q,2) = diviuexact(y,t);
     417             :     }
     418             :   }
     419             :   else
     420             :   { /* rare: r and q left on stack for efficiency */
     421      883606 :     q = cgetg(3,t_FRAC);
     422      883622 :     gel(q,1) = diviiexact(x,r);
     423      883622 :     gel(q,2) = diviiexact(y,r);
     424             :   }
     425     8378070 :   normalize_frac(q); return q;
     426             : }
     427             : 
     428             : /* x t_INT, return x/y in reduced form */
     429             : GEN
     430   133489560 : Qdiviu(GEN x, ulong y)
     431             : {
     432   133489560 :   pari_sp av = avma;
     433             :   ulong r, t;
     434             :   GEN q;
     435             : 
     436   133489560 :   if (y == 1) return icopy(x);
     437   112705104 :   if (!y) pari_err_INV("gdiv",gen_0);
     438   112705104 :   if (equali1(x)) retmkfrac(gen_1, utoipos(y));
     439   105576858 :   q = absdiviu_rem(x,y,&r);
     440   105572657 :   if (!r)
     441             :   {
     442    61871065 :     if (signe(x) < 0) togglesign(q);
     443    61871098 :     return q;
     444             :   }
     445    43701592 :   t = ugcd(y, r); set_avma(av);
     446    43708967 :   if (t == 1) retmkfrac(icopy(x), utoipos(y));
     447    17389729 :   retmkfrac(diviuexact(x,t), utoipos(y / t));
     448             : }
     449             : 
     450             : /* x t_INT, return x/y in reduced form */
     451             : GEN
     452     1561457 : Qdivis(GEN x, long y)
     453             : {
     454     1561457 :   pari_sp av = avma;
     455             :   ulong r, t;
     456             :   long s;
     457             :   GEN q;
     458             : 
     459     1561457 :   if (y > 0) return Qdiviu(x, y);
     460      162430 :   if (!y) pari_err_INV("gdiv",gen_0);
     461      162430 :   s = signe(x);
     462      162430 :   if (!s) return gen_0;
     463      132202 :   if (y < 0) { y = -y; s = -s; }
     464      132202 :   if (y == 1) { x = icopy(x); setsigne(x,s); return x; }
     465      131922 :   if (equali1(x)) retmkfrac(s > 0? gen_1: gen_m1, utoipos(y));
     466      129759 :   q = absdiviu_rem(x,y,&r);
     467      129759 :   if (!r)
     468             :   {
     469       55356 :     if (s < 0) togglesign(q);
     470       55356 :     return q;
     471             :   }
     472       74403 :   t = ugcd(y, r); set_avma(av); q = cgetg(3, t_FRAC);
     473       74403 :   if (t != 1) { x = diviuexact(x,t); y /= t; } else x = icopy(x);
     474       74403 :   gel(q,1) = x; setsigne(x, s);
     475       74403 :   gel(q,2) = utoipos(y); return q;
     476             : }
     477             : 
     478             : /*******************************************************************/
     479             : /*                                                                 */
     480             : /*                          CONJUGATION                            */
     481             : /*                                                                 */
     482             : /*******************************************************************/
     483             : /* lift( conj(Mod(x, y)) ), assuming degpol(y) = 2, degpol(x) < 2 */
     484             : static GEN
     485       18340 : quad_polmod_conj(GEN x, GEN y)
     486             : {
     487             :   GEN z, u, v, a, b;
     488       18340 :   if (typ(x) != t_POL) return gcopy(x);
     489       18340 :   if (varn(x) != varn(y) || degpol(x) <= 0) return RgX_copy(x);
     490       18340 :   a = gel(y,4); u = gel(x,3); /*Mod(ux + v, ax^2 + bx + c)*/
     491       18340 :   b = gel(y,3); v = gel(x,2);
     492       18340 :   z = cgetg(4, t_POL); z[1] = x[1];
     493       18340 :   gel(z,2) = gsub(v, gdiv(gmul(u,b), a));
     494       18340 :   gel(z,3) = gneg(u); return z;
     495             : }
     496             : static GEN
     497       18340 : quad_polmod_norm(GEN x, GEN y)
     498             : {
     499             :   GEN z, u, v, a, b, c;
     500       18340 :   if (typ(x) != t_POL || varn(x) != varn(y) || degpol(x) <= 0) return gsqr(x);
     501       18340 :   a = gel(y,4); u = gel(x,3); /*Mod(ux + v, ax^2 + bx + c)*/
     502       18340 :   b = gel(y,3); v = gel(x,2);
     503       18340 :   c = gel(y,2);
     504       18340 :   z = gmul(u, gsub(gmul(c,u), gmul(b,v)));
     505       18340 :   if (!gequal1(a)) z = gdiv(z, a);
     506       18340 :   return gadd(z, gsqr(v));
     507             : }
     508             : 
     509             : GEN
     510    19580962 : conj_i(GEN x)
     511             : {
     512    19580962 :   switch(typ(x))
     513             :   {
     514     6541357 :     case t_INT: case t_REAL: case t_INTMOD: case t_FRAC: case t_PADIC:
     515     6541357 :       return x;
     516             : 
     517    12875943 :     case t_COMPLEX: return mkcomplex(gel(x,1), gneg(gel(x,2)));
     518             : 
     519         924 :     case t_QUAD:
     520             :     {
     521         924 :       GEN y = cgetg(4,t_QUAD);
     522         924 :       gel(y,1) = gel(x,1);
     523         924 :       gel(y,2) = gequal0(gmael(x,1,3))? gel(x,2)
     524         924 :                                       : gadd(gel(x,2), gel(x,3));
     525         924 :       gel(y,3) = gneg(gel(x,3)); return y;
     526             :     }
     527         350 :     case t_POL: pari_APPLY_pol_normalized(conj_i(gel(x,i)));
     528       31591 :     case t_SER: pari_APPLY_ser_normalized(conj_i(gel(x,i)));
     529             : 
     530      152836 :     case t_RFRAC:
     531             :     case t_VEC:
     532             :     case t_COL:
     533      550625 :     case t_MAT: pari_APPLY_same(conj_i(gel(x,i)));
     534             : 
     535           0 :     case t_POLMOD:
     536             :     {
     537           0 :       GEN X = gel(x,1);
     538           0 :       long d = degpol(X);
     539           0 :       if (d < 2) return x;
     540           0 :       if (d == 2) return mkpolmod(quad_polmod_conj(gel(x,2), X), X);
     541             :     }
     542             :   }
     543           0 :   pari_err_TYPE("gconj",x);
     544             :   return NULL; /* LCOV_EXCL_LINE */
     545             : }
     546             : GEN
     547      133221 : gconj(GEN x)
     548      133221 : { pari_sp av = avma; return gerepilecopy(av, conj_i(x)); }
     549             : 
     550             : GEN
     551          84 : conjvec(GEN x,long prec)
     552             : {
     553             :   long lx, s, i;
     554             :   GEN z;
     555             : 
     556          84 :   switch(typ(x))
     557             :   {
     558           0 :     case t_INT: case t_INTMOD: case t_FRAC:
     559           0 :       return mkcolcopy(x);
     560             : 
     561           0 :     case t_COMPLEX: case t_QUAD:
     562           0 :       z=cgetg(3,t_COL); gel(z,1) = gcopy(x); gel(z,2) = gconj(x); break;
     563             : 
     564          28 :     case t_FFELT:
     565          28 :       return FF_conjvec(x);
     566             : 
     567           0 :     case t_VEC: case t_COL:
     568           0 :       lx = lg(x); z = cgetg(lx,t_MAT);
     569           0 :       if (lx == 1) return z;
     570           0 :       gel(z,1) = conjvec(gel(x,1),prec);
     571           0 :       s = lgcols(z);
     572           0 :       for (i=2; i<lx; i++)
     573             :       {
     574           0 :         gel(z,i) = conjvec(gel(x,i),prec);
     575           0 :         if (lg(gel(z,i)) != s) pari_err_OP("conjvec", gel(z,1), gel(z,i));
     576             :       }
     577           0 :       break;
     578             : 
     579          56 :     case t_POLMOD: {
     580          56 :       GEN T = gel(x,1), r;
     581             :       pari_sp av;
     582             : 
     583          56 :       lx = lg(T);
     584          56 :       if (lx <= 3) return cgetg(1,t_COL);
     585          56 :       x = gel(x,2);
     586         238 :       for (i=2; i<lx; i++)
     587             :       {
     588         189 :         GEN c = gel(T,i);
     589         189 :         switch(typ(c)) {
     590           7 :           case t_INTMOD: {
     591           7 :             GEN p = gel(c,1);
     592             :             pari_sp av;
     593           7 :             if (typ(x) != t_POL) retconst_col(lx-3, Rg_to_Fp(x, p));
     594           7 :             av = avma;
     595           7 :             T = RgX_to_FpX(T,p);
     596           7 :             x = RgX_to_FpX(x, p);
     597           7 :             if (varn(x) != varn(T)) pari_err_VAR("conjvec",x,T);
     598           7 :             z = FpXQC_to_mod(FpXQ_conjvec(x, T , p), T, p);
     599           7 :             return gerepileupto(av, z);
     600             :           }
     601         182 :           case t_INT:
     602         182 :           case t_FRAC: break;
     603           0 :           default: pari_err_TYPE("conjvec [not a rational t_POL]",T);
     604             :         }
     605             :       }
     606          49 :       if (typ(x) != t_POL)
     607             :       {
     608           0 :         if (!is_rational_t(typ(x)))
     609           0 :           pari_err_TYPE("conjvec [not a rational t_POL]",x);
     610           0 :         retconst_col(lx-3, gcopy(x));
     611             :       }
     612          49 :       RgX_check_QX(x,"conjvec");
     613          49 :       av = avma;
     614          49 :       if (varn(x) != varn(T)) pari_err_VAR("conjvec",x,T);
     615          49 :       r = cleanroots(T,prec);
     616          49 :       z = cgetg(lx-2,t_COL);
     617         182 :       for (i=1; i<=lx-3; i++) gel(z,i) = poleval(x, gel(r,i));
     618          49 :       return gerepileupto(av, z);
     619             :     }
     620             : 
     621           0 :     default:
     622           0 :       pari_err_TYPE("conjvec",x);
     623             :       return NULL; /* LCOV_EXCL_LINE */
     624             :   }
     625           0 :   return z;
     626             : }
     627             : 
     628             : /********************************************************************/
     629             : /**                                                                **/
     630             : /**                           ADDITION                             **/
     631             : /**                                                                **/
     632             : /********************************************************************/
     633             : static GEN
     634    24588716 : mkpadic_mod(GEN u, GEN p, GEN pd, long e, long d)
     635    24588716 : { retmkpadic(modii(u, pd), icopy(p), icopy(pd), e, d); }
     636             : 
     637             : /* x, y compatible PADIC, op = add or sub */
     638             : static GEN
     639    19968181 : addsub_pp(GEN x, GEN y, GEN (*op)(GEN,GEN))
     640             : {
     641    19968181 :   pari_sp av = avma;
     642             :   long d, e, r, rx, ry;
     643    19968181 :   GEN u, z, mod, pdx, pdy, ux, uy, p = padic_p(x);
     644             :   int swap;
     645             : 
     646    19968181 :   e = valp(x);
     647    19968181 :   r = valp(y); d = r-e;
     648    19968181 :   if (d < 0) { swap = 1; swap(x,y); e = r; d = -d; } else swap = 0;
     649    19968181 :   pdx = padic_pd(x); ux = padic_u(x);
     650    19968181 :   pdy = padic_pd(y); uy = padic_u(y);
     651    19968181 :   rx = precp(x);
     652    19968181 :   ry = precp(y);
     653    19968181 :   if (d) /* v(x) < v(y) */
     654             :   {
     655    10703429 :     r = d+ry; z = powiu(p,d);
     656    10703492 :     if (r < rx) mod = mulii(z,pdy); else { r = rx; mod = pdx; }
     657    10703510 :     z = mulii(z, uy);
     658    10703479 :     u = swap? op(z, ux): op(ux, z);
     659             :   }
     660             :   else
     661             :   {
     662             :     long c;
     663     9264752 :     if (ry < rx) { r=ry; mod = pdy; } else { r=rx; mod = pdx; }
     664     9264752 :     u = swap? op(uy, ux): op(ux, uy);
     665     9265413 :     if (!signe(u) || (c = Z_pvalrem(u,p,&u)) >= r)
     666             :     {
     667      138198 :       set_avma(av); return zeropadic(p, e+r);
     668             :     }
     669     9127026 :     if (c)
     670             :     {
     671     3351441 :       mod = diviiexact(mod, powiu(p,c));
     672     3351446 :       r -= c;
     673     3351446 :       e += c;
     674             :     }
     675             :   }
     676    19830511 :   return gerepileupto(av, mkpadic_mod(u, p, mod, e, r));
     677             : }
     678             : /* Rg_to_Fp(t_FRAC) without GC */
     679             : static GEN
     680      226838 : Q_to_Fp(GEN x, GEN p)
     681      226838 : { return mulii(gel(x,1), Fp_inv(gel(x,2),p)); }
     682             : /* return x + y, where y t_PADIC and x is a nonzero t_INT or t_FRAC */
     683             : static GEN
     684     4762904 : addQp(GEN x, GEN y)
     685             : {
     686     4762904 :   pari_sp av = avma;
     687     4762904 :   long d, r, e, vy = valp(y), py = precp(y);
     688     4762904 :   GEN q, mod, u, p = padic_p(y);
     689             : 
     690     4762904 :   e = Q_pvalrem(x, p, &x);
     691     4762899 :   d = vy - e; r = d + py;
     692     4762899 :   if (r <= 0) { set_avma(av); return gcopy(y); }
     693     4761015 :   mod = padic_pd(y);
     694     4761015 :   u   = padic_u(y);
     695     4761015 :   if (d > 0)
     696             :   {
     697     1376559 :     q = powiu(p,d);
     698     1376563 :     mod = mulii(mod, q);
     699     1376561 :     if (typ(x) != t_INT) x = Q_to_Fp(x, mod);
     700     1376561 :     u = addii(x,  mulii(u, q));
     701             :   }
     702     3384456 :   else if (d < 0)
     703             :   {
     704      405332 :     q = powiu(p,-d);
     705      405332 :     if (typ(x) != t_INT) x = Q_to_Fp(x, mod);
     706      405332 :     u = addii(u, mulii(x, q));
     707      405332 :     r = py; e = vy;
     708             :   }
     709             :   else
     710             :   {
     711             :     long c;
     712     2979124 :     if (typ(x) != t_INT) x = Q_to_Fp(x, mod);
     713     2979124 :     u = addii(u, x);
     714     2979125 :     if (!signe(u) || (c = Z_pvalrem(u,p,&u)) >= r)
     715             :     {
     716        1204 :       set_avma(av); return zeropadic(p,e+r);
     717             :     }
     718     2977921 :     if (c)
     719             :     {
     720      970645 :       mod = diviiexact(mod, powiu(p,c));
     721      970645 :       r -= c;
     722      970645 :       e += c;
     723             :     }
     724             :   }
     725     4759817 :   return gerepileupto(av, mkpadic_mod(u, p, mod, e, r));
     726             : }
     727             : 
     728             : /* Mod(x,X) + Mod(y,X) */
     729             : #define addsub_polmod_same addsub_polmod_scal
     730             : /* Mod(x,X) +/- Mod(y,Y) */
     731             : static GEN
     732        7203 : addsub_polmod(GEN X, GEN Y, GEN x, GEN y, GEN(*op)(GEN,GEN))
     733             : {
     734        7203 :   long T[3] = { evaltyp(t_POLMOD) | _evallg(3),0,0 };
     735        7203 :   GEN z = cgetg(3,t_POLMOD);
     736        7203 :   long vx = varn(X), vy = varn(Y);
     737        7203 :   if (vx==vy) {
     738             :     pari_sp av;
     739          14 :     gel(z,1) = RgX_gcd(X,Y); av = avma;
     740          14 :     warn_coercion(X,Y,gel(z,1));
     741          14 :     gel(z,2) = gerepileupto(av, gmod(op(x, y), gel(z,1))); return z;
     742             :   }
     743        7189 :   if (varncmp(vx, vy) < 0)
     744        7189 :   { gel(z,1) = RgX_copy(X); gel(T,1) = Y; gel(T,2) = y; y = T; }
     745             :   else
     746           0 :   { gel(z,1) = RgX_copy(Y); gel(T,1) = X; gel(T,2) = x; x = T; }
     747        7189 :   gel(z,2) = op(x, y); return z;
     748             : }
     749             : /* Mod(y, Y) +/- x,  x scalar or polynomial in same var and reduced degree */
     750             : static GEN
     751    13420752 : addsub_polmod_scal(GEN Y, GEN y, GEN x, GEN(*op)(GEN,GEN))
     752    13420752 : { retmkpolmod(degpol(Y)? op(y, x): gen_0, RgX_copy(Y)); }
     753             : 
     754             : /* typ(y) == t_SER, x "scalar" [e.g object in lower variable] */
     755             : static GEN
     756      408094 : add_ser_scal(GEN y, GEN x)
     757             : {
     758             :   long i, v, ly, vy;
     759             :   GEN z;
     760             : 
     761      408094 :   if (isrationalzero(x)) return gcopy(y);
     762      381837 :   ly = lg(y);
     763      381837 :   v = valser(y);
     764      381837 :   if (v < 3-ly) return gcopy(y);
     765             :   /* v + ly >= 3 */
     766      381585 :   if (v < 0)
     767             :   {
     768        1162 :     z = cgetg(ly,t_SER); z[1] = y[1];
     769        3276 :     for (i = 2; i <= 1-v; i++) gel(z,i) = gcopy(gel(y,i));
     770        1162 :     gel(z,i) = gadd(x,gel(y,i)); i++;
     771        3157 :     for (     ; i < ly; i++)   gel(z,i) = gcopy(gel(y,i));
     772        1162 :     return normalizeser(z);
     773             :   }
     774      380423 :   vy = varn(y);
     775      380423 :   if (v > 0)
     776             :   {
     777       19299 :     if (ser_isexactzero(y))
     778        9338 :       return scalarser(ly == 2? x: gadd(x,gel(y,2)), vy, v);
     779        9961 :     y -= v; ly += v;
     780        9961 :     z = cgetg(ly,t_SER);
     781        9961 :     x = gcopy(x);
     782       20559 :     for (i=3; i<=v+1; i++) gel(z,i) = gen_0;
     783             :   }
     784      361124 :   else if (ser_isexactzero(y)) return gcopy(y);
     785             :   else
     786             :   { /* v = 0, ly >= 3 */
     787      361117 :     z = cgetg(ly,t_SER);
     788      361117 :     x = gadd(x, gel(y,2));
     789      361117 :     i = 3;
     790             :   }
     791     1598979 :   for (; i<ly; i++) gel(z,i) = gcopy(gel(y,i));
     792      371078 :   gel(z,2) = x;
     793      371078 :   z[1] = evalsigne(1) | _evalvalser(0) | evalvarn(vy);
     794      371078 :   return gequal0(x)? normalizeser(z): z;
     795             : }
     796             : static long
     797     7161978 : _serprec(GEN x) { return ser_isexactzero(x)? 2: lg(x); }
     798             : /* x,y t_SER in the same variable: x+y */
     799             : static GEN
     800     3581388 : ser_add(GEN x, GEN y)
     801             : {
     802     3581388 :   long i, lx,ly, n = valser(y) - valser(x);
     803             :   GEN z;
     804     3581388 :   if (n < 0) { n = -n; swap(x,y); }
     805             :   /* valser(x) <= valser(y) */
     806     3581388 :   lx = _serprec(x);
     807     3581388 :   if (lx == 2) /* don't lose type information */
     808             :   {
     809         798 :     z = scalarser(gadd(Rg_get_0(x), Rg_get_0(y)), varn(x), 1);
     810         798 :     setvalser(z, valser(x)); return z;
     811             :   }
     812     3580590 :   ly = _serprec(y) + n; if (lx < ly) ly = lx;
     813     3580590 :   if (n)
     814             :   {
     815      107782 :     if (n+2 > lx) return gcopy(x);
     816      107082 :     z = cgetg(ly,t_SER);
     817      819116 :     for (i=2; i<=n+1; i++) gel(z,i) = gcopy(gel(x,i));
     818      506624 :     for (   ; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i-n));
     819             :   } else {
     820     3472808 :     z = cgetg(ly,t_SER);
     821    17589354 :     for (i=2; i < ly; i++) gel(z,i) = gadd(gel(x,i),gel(y,i));
     822             :   }
     823     3579890 :   z[1] = x[1]; return normalizeser(z);
     824             : }
     825             : /* typ(y) == RFRAC, x polynomial in same variable or "scalar" */
     826             : static GEN
     827     8811856 : add_rfrac_scal(GEN y, GEN x)
     828             : {
     829             :   pari_sp av;
     830             :   GEN n;
     831             : 
     832     8811856 :   if (isintzero(x)) return gcopy(y); /* frequent special case */
     833     5138091 :   av = avma; n = gadd(gmul(x, gel(y,2)), gel(y,1));
     834     5138091 :   return gerepileupto(av, gred_rfrac_simple(n, gel(y,2)));
     835             : }
     836             : 
     837             : /* x "scalar", ty != t_MAT and nonscalar */
     838             : static GEN
     839    44613410 : add_scal(GEN y, GEN x, long ty)
     840             : {
     841    44613410 :   switch(ty)
     842             :   {
     843    39701190 :     case t_POL: return RgX_Rg_add(y, x);
     844      408066 :     case t_SER: return add_ser_scal(y, x);
     845     4468778 :     case t_RFRAC: return add_rfrac_scal(y, x);
     846           0 :     case t_COL: return RgC_Rg_add(y, x);
     847       35370 :     case t_VEC:
     848       35370 :       if (isintzero(x)) return gcopy(y);
     849         175 :       break;
     850             :   }
     851         181 :   pari_err_TYPE2("+",x,y);
     852             :   return NULL; /* LCOV_EXCL_LINE */
     853             : }
     854             : 
     855             : /* assumes z = cget(3, t_FRAC) comes first in stack, then a, then b */
     856             : static GEN
     857    14436549 : setfrac(GEN z, GEN a, GEN b)
     858             : {
     859    14436549 :   gel(z,1) = icopy_avma(a, (pari_sp)z);
     860    14436545 :   gel(z,2) = icopy_avma(b, (pari_sp)gel(z,1));
     861    14436629 :   set_avma((pari_sp)gel(z,2)); return z;
     862             : }
     863             : /* z <- a / (b*Q), (Q,a) = 1 */
     864             : static GEN
     865    13628016 : addsub_frac_i(GEN z, GEN Q, GEN a, GEN b)
     866             : {
     867    13628016 :   GEN q = Qdivii(a, b); /* != 0 */
     868    13628068 :   if (typ(q) == t_INT)
     869             :   {
     870     1951427 :     gel(z,1) = gerepileuptoint((pari_sp)Q, q);
     871     1951427 :     gel(z,2) = Q; return z;
     872             :   }
     873    11676641 :   return setfrac(z, gel(q,1), mulii(gel(q,2), Q));
     874             : }
     875             : static GEN
     876    27392414 : addsub_frac(GEN x, GEN y, GEN (*op)(GEN,GEN))
     877             : {
     878    27392414 :   GEN x1 = gel(x,1), x2 = gel(x,2);
     879    27392414 :   GEN y1 = gel(y,1), y2 = gel(y,2), z, Q, q, r, n, delta;
     880    27392414 :   int s = cmpii(x2, y2);
     881             : 
     882             :   /* common denominator: (x1 op y1) / x2 */
     883    27392421 :   if (!s)
     884             :   {
     885     9465052 :     pari_sp av = avma;
     886     9465052 :     return gerepileupto(av, Qdivii(op(x1, y1), x2));
     887             :   }
     888    17927369 :   z = cgetg(3, t_FRAC);
     889    17927376 :   if (s < 0)
     890             :   {
     891    10262463 :     Q = dvmdii(y2, x2, &r);
     892             :     /* y2 = Q x2: 1/x2 . (Q x1 op y1)/Q, where latter is in coprime form */
     893    10262418 :     if (r == gen_0) return addsub_frac_i(z, Q, op(mulii(Q,x1), y1), x2);
     894     2740080 :     delta = gcdii(x2,r);
     895             :   }
     896             :   else
     897             :   {
     898     7664913 :     Q = dvmdii(x2, y2, &r);
     899             :     /* x2 = Q y2: 1/y2 . (x1 op Q y1)/Q, where latter is in coprime form */
     900     7664919 :     if (r == gen_0) return addsub_frac_i(z, Q, op(x1, mulii(Q,y1)), y2);
     901     1559185 :     delta = gcdii(y2,r);
     902             :   }
     903             :   /* delta = gcd(x2,y2) */
     904     4299304 :   if (equali1(delta))
     905             :   { /* numerator is nonzero */
     906     1539323 :     gel(z,1) = gerepileuptoint((pari_sp)z, op(mulii(x1,y2), mulii(y1,x2)));
     907     1539323 :     gel(z,2) = mulii(x2,y2); return z;
     908             :   }
     909     2759979 :   x2 = diviiexact(x2,delta);
     910     2759977 :   y2 = diviiexact(y2,delta);
     911     2759973 :   n = op(mulii(x1,y2), mulii(y1,x2)); /* != 0 */
     912     2759980 :   q = dvmdii(n, delta, &r);
     913     2759975 :   if (r == gen_0) return setfrac(z, q, mulii(x2, y2));
     914     2506735 :   r = gcdii(delta, r);
     915     2506739 :   if (!equali1(r)) { n = diviiexact(n, r); delta = diviiexact(delta, r); }
     916     2506739 :   return setfrac(z, n, mulii(mulii(x2, y2), delta));
     917             : }
     918             : 
     919             : /* assume x2, y2 are t_POLs in the same variable */
     920             : static GEN
     921     3018187 : add_rfrac(GEN x, GEN y)
     922             : {
     923     3018187 :   pari_sp av = avma;
     924     3018187 :   GEN x1 = gel(x,1), x2 = gel(x,2);
     925     3018187 :   GEN y1 = gel(y,1), y2 = gel(y,2), q, r, n, d, delta;
     926             : 
     927     3018187 :   delta = RgX_gcd(x2,y2);
     928     3018187 :   if (!degpol(delta))
     929             :   {
     930         644 :     n = simplify_shallow( gadd(gmul(x1,y2), gmul(y1,x2)) );
     931         644 :     d = RgX_mul(x2, y2);
     932         644 :     return gerepileupto(av, gred_rfrac_simple(n, d));
     933             :   }
     934     3017543 :   x2 = RgX_div(x2,delta);
     935     3017543 :   y2 = RgX_div(y2,delta);
     936     3017543 :   n = gadd(gmul(x1,y2), gmul(y1,x2));
     937     3017543 :   if (!signe(n))
     938             :   {
     939      721418 :     n = simplify_shallow(n);
     940      721418 :     if (isexactzero(n))
     941             :     {
     942      721411 :       if (isrationalzero(n)) { set_avma(av); return zeropol(varn(x2)); }
     943           7 :       return gerepileupto(av, scalarpol(n, varn(x2)));
     944             :     }
     945           7 :     return gerepilecopy(av, mkrfrac(n, RgX_mul(gel(x,2),y2)));
     946             :   }
     947     2296125 :   if (degpol(n) == 0)
     948     1150575 :     return gerepileupto(av, gred_rfrac_simple(gel(n,2), RgX_mul(gel(x,2),y2)));
     949     1145550 :   q = RgX_divrem(n, delta, &r); /* we want gcd(n,delta) */
     950     1145550 :   if (isexactzero(r))
     951             :   {
     952             :     GEN z;
     953      227982 :     d = RgX_mul(x2, y2);
     954             :     /* "constant" denominator ? */
     955      227982 :     z = lg(d) == 3? RgX_Rg_div(q, gel(d,2)): gred_rfrac_simple(q, d);
     956      227982 :     return gerepileupto(av, z);
     957             :   }
     958      917568 :   r = RgX_gcd(delta, r);
     959      917568 :   if (degpol(r))
     960             :   {
     961      160642 :     n = RgX_div(n, r);
     962      160642 :     d = RgX_mul(RgX_mul(x2,y2), RgX_div(delta, r));
     963             :   }
     964             :   else
     965      756926 :     d = RgX_mul(gel(x,2), y2);
     966      917568 :   return gerepileupto(av, gred_rfrac_simple(n, d));
     967             : }
     968             : 
     969             : GEN
     970  6054199207 : gadd(GEN x, GEN y)
     971             : {
     972  6054199207 :   long tx = typ(x), ty = typ(y), vx, vy, lx, i, l;
     973             :   pari_sp av, tetpil;
     974             :   GEN z, p1;
     975             : 
     976  6054199207 :   if (tx == ty) switch(tx) /* shortcut to generic case */
     977             :   {
     978  2996110528 :     case t_INT: return addii(x,y);
     979  1639294735 :     case t_REAL: return addrr(x,y);
     980     1381498 :     case t_INTMOD:  { GEN X = gel(x,1), Y = gel(y,1);
     981     1381498 :       z = cgetg(3,t_INTMOD);
     982     1381498 :       if (X==Y || equalii(X,Y))
     983     1381484 :         return add_intmod_same(z, X, gel(x,2), gel(y,2));
     984          14 :       gel(z,1) = gcdii(X,Y);
     985          14 :       warn_coercion(X,Y,gel(z,1));
     986          14 :       av = avma; p1 = addii(gel(x,2),gel(y,2));
     987          14 :       gel(z,2) = gerepileuptoint(av, remii(p1, gel(z,1))); return z;
     988             :     }
     989    22322991 :     case t_FRAC: return addsub_frac(x,y,addii);
     990   354930823 :     case t_COMPLEX: z = cgetg(3,t_COMPLEX);
     991   354447787 :       gel(z,2) = gadd(gel(x,2),gel(y,2));
     992   355214738 :       if (isintzero(gel(z,2)))
     993             :       {
     994      521178 :         set_avma((pari_sp)(z+3));
     995      521178 :         return gadd(gel(x,1),gel(y,1));
     996             :       }
     997   354577510 :       gel(z,1) = gadd(gel(x,1),gel(y,1));
     998   354570480 :       return z;
     999    14530083 :     case t_PADIC:
    1000    14530083 :       if (!equalii(padic_p(x), padic_p(y))) pari_err_OP("+",x,y);
    1001    14529844 :       return addsub_pp(x,y, addii);
    1002         672 :     case t_QUAD: z = cgetg(4,t_QUAD);
    1003         672 :       if (!ZX_equal(gel(x,1),gel(y,1))) pari_err_OP("+",x,y);
    1004         672 :       gel(z,1) = ZX_copy(gel(x,1));
    1005         672 :       gel(z,2) = gadd(gel(x,2),gel(y,2));
    1006         672 :       gel(z,3) = gadd(gel(x,3),gel(y,3)); return z;
    1007    10156362 :     case t_POLMOD:
    1008    10156362 :       if (RgX_equal_var(gel(x,1), gel(y,1)))
    1009    10149194 :         return addsub_polmod_same(gel(x,1), gel(x,2), gel(y,2), &gadd);
    1010        7168 :       return addsub_polmod(gel(x,1), gel(y,1), gel(x,2), gel(y,2), &gadd);
    1011    13133535 :     case t_FFELT: return FF_add(x,y);
    1012   139553342 :     case t_POL:
    1013   139553342 :       vx = varn(x);
    1014   139553342 :       vy = varn(y);
    1015   139553342 :       if (vx != vy) {
    1016     4856925 :         if (varncmp(vx, vy) < 0) return RgX_Rg_add(x, y);
    1017     1798442 :         else                     return RgX_Rg_add(y, x);
    1018             :       }
    1019   134696417 :       return RgX_add(x, y);
    1020     3576033 :     case t_SER:
    1021     3576033 :       vx = varn(x);
    1022     3576033 :       vy = varn(y);
    1023     3576033 :       if (vx != vy) {
    1024          28 :         if (varncmp(vx, vy) < 0) return add_ser_scal(x, y);
    1025          21 :         else                     return add_ser_scal(y, x);
    1026             :       }
    1027     3576005 :       return ser_add(x, y);
    1028     4310088 :     case t_RFRAC:
    1029     4310088 :       vx = varn(gel(x,2));
    1030     4310088 :       vy = varn(gel(y,2));
    1031     4310088 :       if (vx != vy) {
    1032     1291901 :         if (varncmp(vx, vy) < 0) return add_rfrac_scal(x, y);
    1033      538397 :         else                     return add_rfrac_scal(y, x);
    1034             :       }
    1035     3018187 :       return add_rfrac(x,y);
    1036     2060586 :     case t_VEC:
    1037     2060586 :       if (lg(y) != lg(x)) pari_err_OP("+",x,y);
    1038     2060586 :       return RgV_add(x,y);
    1039     1103269 :     case t_COL:
    1040     1103269 :       if (lg(y) != lg(x)) pari_err_OP("+",x,y);
    1041     1103269 :       return RgC_add(x,y);
    1042      671412 :     case t_MAT:
    1043      671412 :       lx = lg(x);
    1044      671412 :       if (lg(y) != lx) pari_err_OP("+",x,y);
    1045      671412 :       if (lx == 1) return cgetg(1, t_MAT);
    1046      671412 :       if (lgcols(y) != lgcols(x)) pari_err_OP("+",x,y);
    1047      671405 :       return RgM_add(x,y);
    1048             : 
    1049           0 :     default: pari_err_TYPE2("+",x,y);
    1050             :   }
    1051             :   /* tx != ty */
    1052   855508699 :   if (tx > ty) { swap(x,y); lswap(tx,ty); }
    1053             : 
    1054   855508699 :   if (is_const_t(ty)) switch(tx) /* tx < ty, is_const_t(tx) && is_const_t(ty) */
    1055             :   {
    1056   740197503 :     case t_INT:
    1057             :       switch(ty)
    1058             :       {
    1059   436865370 :         case t_REAL: return addir(x,y);
    1060     2139125 :         case t_INTMOD:
    1061     2139125 :           z = cgetg(3, t_INTMOD);
    1062     2139125 :           return add_intmod_same(z, gel(y,1), gel(y,2), modii(x, gel(y,1)));
    1063    44124282 :         case t_FRAC: z = cgetg(3,t_FRAC);
    1064    44124274 :           gel(z,1) = gerepileuptoint((pari_sp)z, addii(gel(y,1), mulii(gel(y,2),x)));
    1065    44124239 :           gel(z,2) = icopy(gel(y,2)); return z;
    1066   250148121 :         case t_COMPLEX: return addRc(x, y);
    1067     5597444 :         case t_PADIC:
    1068     5597444 :           if (!signe(x)) return gcopy(y);
    1069     4533166 :           return addQp(x,y);
    1070        1113 :         case t_QUAD: return addRq(x, y);
    1071     1365080 :         case t_FFELT: return FF_Z_add(y,x);
    1072             :       }
    1073             : 
    1074             :     case t_REAL:
    1075    57123607 :       switch(ty)
    1076             :       {
    1077    13536632 :         case t_FRAC:
    1078    13536632 :           if (!signe(gel(y,1))) return rcopy(x);
    1079    13536632 :           if (!signe(x))
    1080             :           {
    1081       12084 :             lx = expi(gel(y,1)) - expi(gel(y,2)) - expo(x);
    1082       12084 :             return lx <= 0? rcopy(x): fractor(y, nbits2prec(lx));
    1083             :           }
    1084    13524548 :           av=avma; z=addir(gel(y,1),mulir(gel(y,2),x)); tetpil=avma;
    1085    13523169 :           return gerepile(av,tetpil,divri(z,gel(y,2)));
    1086    43586905 :         case t_COMPLEX: return addRc(x, y);
    1087          70 :         case t_QUAD: return gequal0(y)? rcopy(x): addqf(y, x, realprec(x));
    1088             : 
    1089           0 :         default: pari_err_TYPE2("+",x,y);
    1090             :       }
    1091             : 
    1092       17640 :     case t_INTMOD:
    1093             :       switch(ty)
    1094             :       {
    1095       17507 :         case t_FRAC: { GEN X = gel(x,1);
    1096       17507 :           z = cgetg(3, t_INTMOD);
    1097       17507 :           p1 = Fp_div(gel(y,1), gel(y,2), X);
    1098       17507 :           return add_intmod_same(z, X, p1, gel(x,2));
    1099             :         }
    1100          14 :         case t_FFELT:
    1101          14 :           if (!equalii(gel(x,1),FF_p_i(y)))
    1102           0 :             pari_err_OP("+",x,y);
    1103          14 :           return FF_Z_add(y,gel(x,2));
    1104          91 :         case t_COMPLEX: return addRc(x, y);
    1105           0 :         case t_PADIC: { GEN X = gel(x,1);
    1106           0 :           z = cgetg(3, t_INTMOD);
    1107           0 :           return add_intmod_same(z, X, gel(x,2), padic_to_Fp(y, X));
    1108             :         }
    1109          28 :         case t_QUAD: return addRq(x, y);
    1110             :       }
    1111             : 
    1112             :     case t_FRAC:
    1113             :       switch (ty)
    1114             :       {
    1115    10299335 :         case t_COMPLEX: return addRc(x, y);
    1116      229736 :         case t_PADIC:
    1117      229736 :           if (!signe(gel(x,1))) return gcopy(y);
    1118      229736 :           return addQp(x,y);
    1119         133 :         case t_QUAD: return addRq(x, y);
    1120        1337 :         case t_FFELT: return FF_Q_add(y, x);
    1121             :       }
    1122             : 
    1123             :     case t_FFELT:
    1124           0 :       pari_err_TYPE2("+",x,y);
    1125             : 
    1126          35 :     case t_COMPLEX:
    1127             :       switch(ty)
    1128             :       {
    1129          28 :         case t_PADIC:
    1130          28 :           return Zp_nosquare_m1(padic_p(y))? addRc(y, x): addTp(x, y);
    1131           7 :         case t_QUAD:
    1132           7 :           lx = precision(x); if (!lx) pari_err_OP("+",x,y);
    1133           7 :           return gequal0(y)? gcopy(x): addqf(y, x, lx);
    1134             :       }
    1135             : 
    1136             :     case t_PADIC: /* ty == t_QUAD */
    1137           0 :       return (kro_quad(y, padic_p(x)) == -1)? addRq(x, y): addTp(y, x);
    1138             :   }
    1139             :   /* tx < ty, !is_const_t(y) */
    1140    50284787 :   switch(ty)
    1141             :   {
    1142        7693 :     case t_MAT:
    1143        7693 :       if (is_matvec_t(tx)) pari_err_TYPE2("+",x,y);
    1144        7693 :       if (isrationalzero(x)) return gcopy(y);
    1145        7616 :       return RgM_Rg_add(y, x);
    1146      185446 :     case t_COL:
    1147      185446 :       if (tx == t_VEC) pari_err_TYPE2("+",x,y);
    1148      185446 :       return RgC_Rg_add(y, x);
    1149     2391609 :     case t_POLMOD: /* is_const_t(tx) in this case */
    1150     2391609 :       return addsub_polmod_scal(gel(y,1), gel(y,2), x, &gadd);
    1151             :   }
    1152    47700039 :   if (is_scalar_t(tx))  {
    1153    44620202 :     if (tx == t_POLMOD)
    1154             :     {
    1155      109851 :       vx = varn(gel(x,1));
    1156      109851 :       vy = gvar(y);
    1157      109851 :       if (vx == vy) y = gmod(y, gel(x,1)); /* error if ty == t_SER */
    1158             :       else
    1159       80850 :         if (varncmp(vx,vy) > 0) return add_scal(y, x, ty);
    1160       29540 :       return addsub_polmod_scal(gel(x,1), gel(x,2), y, &gadd);
    1161             :     }
    1162    44510351 :     return add_scal(y, x, ty);
    1163             :   }
    1164             :   /* x and y are not scalars, ty != t_MAT */
    1165     3079839 :   vx = gvar(x);
    1166     3079844 :   vy = gvar(y);
    1167     3079844 :   if (vx != vy) { /* x or y is treated as a scalar */
    1168       22759 :     if (is_vec_t(tx) || is_vec_t(ty)) pari_err_TYPE2("+",x,y);
    1169       32417 :     return (varncmp(vx, vy) < 0)? add_scal(x, y, tx)
    1170       32417 :                                 : add_scal(y, x, ty);
    1171             :   }
    1172             :   /* vx = vy */
    1173     3057085 :   switch(tx)
    1174             :   {
    1175     3056588 :     case t_POL:
    1176             :       switch (ty)
    1177             :       {
    1178        5411 :         case t_SER:
    1179        5411 :           if (lg(x) == 2) return gcopy(y);
    1180        5390 :           i = RgX_val(x); if (i == LONG_MAX) i = 0; /* e.g. x = Mod(0,3)*x^0 */
    1181        5390 :           i = lg(y) + valser(y) - i;
    1182        5390 :           if (i < 3) return gcopy(y);
    1183        5383 :           p1 = RgX_to_ser(x,i); y = ser_add(p1,y);
    1184        5383 :           settyp(p1, t_VECSMALL); /* p1 left on stack */
    1185        5383 :           return y;
    1186             : 
    1187     3051177 :         case t_RFRAC: return add_rfrac_scal(y, x);
    1188             :       }
    1189           0 :       break;
    1190             : 
    1191         497 :     case t_SER:
    1192         497 :       if (ty == t_RFRAC)
    1193             :       {
    1194             :         long vn, vd;
    1195         497 :         av = avma;
    1196         497 :         vn = gval(gel(y,1), vy);
    1197         497 :         vd = RgX_valrem_inexact(gel(y,2), NULL);
    1198         497 :         if (vd == LONG_MAX) pari_err_INV("gadd", gel(y,2));
    1199             : 
    1200         497 :         l = lg(x) + valser(x) - (vn - vd);
    1201         497 :         if (l < 3) { set_avma(av); return gcopy(x); }
    1202         497 :         return gerepileupto(av, gadd(x, rfrac_to_ser_i(y, l)));
    1203             :       }
    1204           0 :       break;
    1205             :   }
    1206           0 :   pari_err_TYPE2("+",x,y);
    1207             :   return NULL; /* LCOV_EXCL_LINE */
    1208             : }
    1209             : 
    1210             : GEN
    1211   270572700 : gaddsg(long x, GEN y)
    1212             : {
    1213   270572700 :   long ty = typ(y);
    1214             :   GEN z;
    1215             : 
    1216   270572700 :   switch(ty)
    1217             :   {
    1218   124842580 :     case t_INT:  return addsi(x,y);
    1219   119175659 :     case t_REAL: return addsr(x,y);
    1220       12003 :     case t_INTMOD:
    1221       12003 :       z = cgetg(3, t_INTMOD);
    1222       12003 :       return add_intmod_same(z, gel(y,1), gel(y,2), modsi(x, gel(y,1)));
    1223    15247831 :     case t_FRAC: z = cgetg(3,t_FRAC);
    1224    15247831 :       gel(z,1) = gerepileuptoint((pari_sp)z, addii(gel(y,1), mulis(gel(y,2),x)));
    1225    15247831 :       gel(z,2) = icopy(gel(y,2)); return z;
    1226     8164429 :     case t_COMPLEX:
    1227     8164429 :       z = cgetg(3, t_COMPLEX);
    1228     8164429 :       gel(z,1) = gaddsg(x, gel(y,1));
    1229     8164429 :       gel(z,2) = gcopy(gel(y,2)); return z;
    1230             : 
    1231     3130198 :     default: return gadd(stoi(x), y);
    1232             :   }
    1233             : }
    1234             : 
    1235             : GEN
    1236     3167908 : gsubsg(long x, GEN y)
    1237             : {
    1238             :   GEN z, a, b;
    1239             :   pari_sp av;
    1240             : 
    1241     3167908 :   switch(typ(y))
    1242             :   {
    1243      276674 :     case t_INT:  return subsi(x,y);
    1244     1225193 :     case t_REAL: return subsr(x,y);
    1245          56 :     case t_INTMOD:
    1246          56 :       z = cgetg(3, t_INTMOD); a = gel(y,1); b = gel(y,2);
    1247          56 :       return add_intmod_same(z, a, Fp_neg(b,a), modsi(x, a));
    1248      732332 :     case t_FRAC: z = cgetg(3,t_FRAC); a = gel(y,1); b = gel(y,2);
    1249      732332 :       gel(z,1) = gerepileuptoint((pari_sp)z, subii(mulis(b,x), a));
    1250      732332 :       gel(z,2) = icopy(gel(y,2)); return z;
    1251      896539 :     case t_COMPLEX:
    1252      896539 :       z = cgetg(3, t_COMPLEX);
    1253      896539 :       gel(z,1) = gsubsg(x, gel(y,1));
    1254      896539 :       gel(z,2) = gneg(gel(y,2)); return z;
    1255             :   }
    1256       37114 :   av = avma;
    1257       37114 :   return gerepileupto(av, gadd(stoi(x), gneg_i(y)));
    1258             : }
    1259             : 
    1260             : /********************************************************************/
    1261             : /**                                                                **/
    1262             : /**                          SUBTRACTION                           **/
    1263             : /**                                                                **/
    1264             : /********************************************************************/
    1265             : 
    1266             : GEN
    1267  2847717119 : gsub(GEN x, GEN y)
    1268             : {
    1269  2847717119 :   long tx = typ(x), ty = typ(y);
    1270             :   pari_sp av;
    1271             :   GEN z;
    1272  2847717119 :   if (tx == ty) switch(tx) /* shortcut to generic case */
    1273             :   {
    1274  2068756742 :     case t_INT: return subii(x,y);
    1275   577711543 :     case t_REAL: return subrr(x,y);
    1276     1053735 :     case t_INTMOD:  { GEN p1, X = gel(x,1), Y = gel(y,1);
    1277     1053735 :       z = cgetg(3,t_INTMOD);
    1278     1053735 :       if (X==Y || equalii(X,Y))
    1279     1053721 :         return sub_intmod_same(z, X, gel(x,2), gel(y,2));
    1280          14 :       gel(z,1) = gcdii(X,Y);
    1281          14 :       warn_coercion(X,Y,gel(z,1));
    1282          14 :       av = avma; p1 = subii(gel(x,2),gel(y,2));
    1283          14 :       gel(z,2) = gerepileuptoint(av, modii(p1, gel(z,1))); return z;
    1284             :     }
    1285     5069436 :     case t_FRAC: return addsub_frac(x,y, subii);
    1286   102061248 :     case t_COMPLEX: z = cgetg(3,t_COMPLEX);
    1287   102098647 :       gel(z,2) = gsub(gel(x,2),gel(y,2));
    1288   102003798 :       if (isintzero(gel(z,2)))
    1289             :       {
    1290       21301 :         set_avma((pari_sp)(z+3));
    1291       21301 :         return gsub(gel(x,1),gel(y,1));
    1292             :       }
    1293   101987453 :       gel(z,1) = gsub(gel(x,1),gel(y,1));
    1294   101968304 :       return z;
    1295     5438381 :     case t_PADIC:
    1296     5438381 :       if (!equalii(padic_p(x), padic_p(y))) pari_err_OP("+",x,y);
    1297     5438381 :       return addsub_pp(x,y, subii);
    1298         168 :     case t_QUAD: z = cgetg(4,t_QUAD);
    1299         168 :       if (!ZX_equal(gel(x,1),gel(y,1))) pari_err_OP("+",x,y);
    1300         168 :       gel(z,1) = ZX_copy(gel(x,1));
    1301         168 :       gel(z,2) = gsub(gel(x,2),gel(y,2));
    1302         168 :       gel(z,3) = gsub(gel(x,3),gel(y,3)); return z;
    1303      850444 :     case t_POLMOD:
    1304      850444 :       if (RgX_equal_var(gel(x,1), gel(y,1)))
    1305      850409 :         return addsub_polmod_same(gel(x,1), gel(x,2), gel(y,2), &gsub);
    1306          35 :       return addsub_polmod(gel(x,1), gel(y,1), gel(x,2), gel(y,2), &gsub);
    1307      240371 :     case t_FFELT: return FF_sub(x,y);
    1308    18420653 :     case t_POL: {
    1309    18420653 :       long vx = varn(x);
    1310    18420653 :       long vy = varn(y);
    1311    18420653 :       if (vx != vy) {
    1312     1146765 :         if (varncmp(vx, vy) < 0) return RgX_Rg_sub(x, y);
    1313       84813 :         else                     return Rg_RgX_sub(x, y);
    1314             :       }
    1315    17273888 :       return RgX_sub(x, y);
    1316             :     }
    1317      299497 :     case t_VEC:
    1318      299497 :       if (lg(y) != lg(x)) pari_err_OP("+",x,y);
    1319      299497 :       return RgV_sub(x,y);
    1320     3569444 :     case t_COL:
    1321     3569444 :       if (lg(y) != lg(x)) pari_err_OP("+",x,y);
    1322     3569444 :       return RgC_sub(x,y);
    1323      252495 :     case t_MAT: {
    1324      252495 :       long lx = lg(x);
    1325      252495 :       if (lg(y) != lx) pari_err_OP("+",x,y);
    1326      252497 :       if (lx == 1) return cgetg(1, t_MAT);
    1327      172867 :       if (lgcols(y) != lgcols(x)) pari_err_OP("+",x,y);
    1328      172868 :       return RgM_sub(x,y);
    1329             :     }
    1330     2448296 :     case t_RFRAC: case t_SER: break;
    1331             : 
    1332           0 :     default: pari_err_TYPE2("+",x,y);
    1333             :   }
    1334    65644966 :   av = avma;
    1335    65644966 :   return gerepileupto(av, gadd(x,gneg_i(y)));
    1336             : }
    1337             : 
    1338             : /********************************************************************/
    1339             : /**                                                                **/
    1340             : /**                        MULTIPLICATION                          **/
    1341             : /**                                                                **/
    1342             : /********************************************************************/
    1343             : static GEN
    1344      316204 : mul_ser_scal(GEN x, GEN t)
    1345             : {
    1346      316204 :   if (isexactzero(t)) return gmul(Rg_get_0(x), t);
    1347      312956 :   if (isint1(t)) return gcopy(x);
    1348      262388 :   if (ser_isexactzero(x))
    1349             :   {
    1350         378 :     GEN z = scalarser(lg(x) == 2? Rg_get_0(t): gmul(gel(x,2), t), varn(x), 1);
    1351         378 :     setvalser(z, valser(x)); return z;
    1352             :   }
    1353     3442418 :   pari_APPLY_ser(gmul(gel(x,i), t));
    1354             : }
    1355             : /* (n/d) * x, x "scalar" or polynomial in the same variable as d
    1356             :  * [n/d a valid RFRAC]  */
    1357             : static GEN
    1358    10416649 : mul_rfrac_scal(GEN n, GEN d, GEN x)
    1359             : {
    1360    10416649 :   pari_sp av = avma;
    1361             :   GEN z;
    1362             : 
    1363    10416649 :   switch(typ(x))
    1364             :   {
    1365          21 :     case t_PADIC:
    1366          21 :       n = gmul(n, x);
    1367          21 :       d = gcvtop(d, padic_p(x), signe(padic_u(x))? precp(x): 1);
    1368          21 :       return gerepileupto(av, gdiv(n,d));
    1369             : 
    1370         623 :     case t_INTMOD: case t_POLMOD:
    1371         623 :       n = gmul(n, x);
    1372         623 :       d = gmul(d, gmodulo(gen_1, gel(x,1)));
    1373         623 :       return gerepileupto(av, gdiv(n,d));
    1374             :   }
    1375    10416005 :   z = gred_rfrac2(x, d);
    1376    10416005 :   n = simplify_shallow(n);
    1377    10416005 :   if (typ(z) == t_RFRAC)
    1378             :   {
    1379     7897143 :     n = gmul(gel(z,1), n);
    1380     7897143 :     d = gel(z,2);
    1381     7897143 :     if (typ(n) == t_POL && varncmp(varn(n), varn(d)) < 0)
    1382           0 :       z = RgX_Rg_div(n, d);
    1383             :     else
    1384     7897143 :       z = gred_rfrac_simple(n, d);
    1385             :   }
    1386             :   else
    1387     2518862 :     z = gmul(z, n);
    1388    10416005 :   return gerepileupto(av, z);
    1389             : }
    1390             : static GEN
    1391   150282286 : mul_scal(GEN y, GEN x, long ty)
    1392             : {
    1393   150282286 :   switch(ty)
    1394             :   {
    1395   141387896 :     case t_POL:
    1396   141387896 :       if (lg(y) == 2) return scalarpol(gmul(gen_0,x), varn(y));
    1397   129469243 :       return RgX_Rg_mul(y, x);
    1398      308063 :     case t_SER: return mul_ser_scal(y, x);
    1399     8586319 :     case t_RFRAC: return mul_rfrac_scal(gel(y,1),gel(y,2), x);
    1400          14 :     case t_QFB:
    1401          14 :       if (typ(x) == t_INT && gequal1(x)) return gcopy(y); /* fall through */
    1402             :   }
    1403          12 :   pari_err_TYPE2("*",x,y);
    1404             :   return NULL; /* LCOV_EXCL_LINE */
    1405             : }
    1406             : static GEN
    1407     7461465 : mul_self_scal(GEN x, GEN y)
    1408   633932890 : { pari_APPLY_same(gmul(y,gel(x,i))); }
    1409             : 
    1410             : static GEN
    1411      160872 : mul_gen_rfrac(GEN X, GEN Y)
    1412             : {
    1413      160872 :   GEN y1 = gel(Y,1), y2 = gel(Y,2);
    1414      160872 :   long vx = gvar(X), vy = varn(y2);
    1415      166633 :   return (varncmp(vx, vy) <= 0)? mul_scal(Y, X, typ(Y)):
    1416        5761 :                                  gred_rfrac_simple(gmul(y1,X), y2);
    1417             : }
    1418             : /* (x1/x2) * (y1/y2) */
    1419             : static GEN
    1420     7907688 : mul_rfrac(GEN x1, GEN x2, GEN y1, GEN y2)
    1421             : {
    1422             :   GEN z, X, Y;
    1423     7907688 :   pari_sp av = avma;
    1424             : 
    1425     7907688 :   X = gred_rfrac2(x1, y2);
    1426     7907688 :   Y = gred_rfrac2(y1, x2);
    1427     7907688 :   if (typ(X) == t_RFRAC)
    1428             :   {
    1429     6628536 :     if (typ(Y) == t_RFRAC) {
    1430     6562939 :       x1 = gel(X,1);
    1431     6562939 :       x2 = gel(X,2);
    1432     6562939 :       y1 = gel(Y,1);
    1433     6562939 :       y2 = gel(Y,2);
    1434     6562939 :       z = gred_rfrac_simple(gmul(x1,y1), gmul(x2,y2));
    1435             :     } else
    1436       65597 :       z = mul_gen_rfrac(Y, X);
    1437             :   }
    1438     1279152 :   else if (typ(Y) == t_RFRAC)
    1439       95275 :     z = mul_gen_rfrac(X, Y);
    1440             :   else
    1441     1183877 :     z = gmul(X, Y);
    1442     7907688 :   return gerepileupto(av, z);
    1443             : }
    1444             : /* (x1/x2) /y2, x2 and y2 are t_POL in the same variable */
    1445             : static GEN
    1446      262512 : div_rfrac_pol(GEN x1, GEN x2, GEN y2)
    1447             : {
    1448      262512 :   pari_sp av = avma;
    1449      262512 :   GEN X = gred_rfrac2(x1, y2);
    1450      262512 :   if (typ(X) == t_RFRAC && varn(gel(X,2)) == varn(x2))
    1451             :   {
    1452      255932 :     x2 = RgX_mul(gel(X,2), x2);
    1453      255932 :     x1 = gel(X,1);
    1454             :   }
    1455             :   else
    1456        6580 :     x1 = X;
    1457      262512 :   return gerepileupto(av, gred_rfrac_simple(x1, x2));
    1458             : }
    1459             : 
    1460             : /* Mod(y, Y) * x,  assuming x scalar */
    1461             : static GEN
    1462     3428803 : mul_polmod_scal(GEN Y, GEN y, GEN x)
    1463     3428803 : { retmkpolmod(gmul(x,y), RgX_copy(Y)); }
    1464             : 
    1465             : /* cf mulqq */
    1466             : static GEN
    1467     5860776 : quad_polmod_mul(GEN P, GEN x, GEN y)
    1468             : {
    1469     5860776 :   GEN T = cgetg(4, t_POL), b = gel(P,3), c = gel(P,2), p1, p2, p3, p4;
    1470     5860776 :   pari_sp tetpil, av = avma;
    1471     5860776 :   T[1] = x[1];
    1472     5860776 :   p2 = gmul(gel(x,2), gel(y,2));
    1473     5860776 :   p3 = gmul(gel(x,3), gel(y,3));
    1474     5860776 :   p1 = gmul(gneg_i(c),p3);
    1475             :   /* operands are usually small: gadd ~ gmul and Karatsuba is a waste */
    1476     5860776 :   if (typ(b) == t_INT)
    1477             :   {
    1478     5860755 :     if (signe(b))
    1479             :     {
    1480     4284989 :       p4 = gadd(gmul(gel(x,2), gel(y,3)), gmul(gel(x,3), gel(y,2)));
    1481     4284989 :       if (is_pm1(b))
    1482             :       {
    1483     4284331 :         if (signe(b) > 0) p3 = gneg(p3);
    1484             :       }
    1485             :       else
    1486         658 :         p3 = gmul(negi(b), p3);
    1487             :     }
    1488             :     else
    1489             :     {
    1490     1575766 :       p3 = gmul(gel(x,2),gel(y,3));
    1491     1575766 :       p4 = gmul(gel(x,3),gel(y,2));
    1492             :     }
    1493             :   }
    1494             :   else
    1495             :   {
    1496          21 :     p4 = gadd(gmul(gel(x,2), gel(y,3)), gmul(gel(x,3), gel(y,2)));
    1497          21 :     p3 = gmul(gneg_i(b), p3);
    1498             :   }
    1499     5860776 :   tetpil = avma;
    1500     5860776 :   gel(T,2) = gadd(p2, p1);
    1501     5860776 :   gel(T,3) = gadd(p4, p3);
    1502     5860776 :   gerepilecoeffssp(av,tetpil,T+2,2);
    1503     5860776 :   return normalizepol_lg(T,4);
    1504             : }
    1505             : /* Mod(x,T) * Mod(y,T) */
    1506             : static GEN
    1507     8565277 : mul_polmod_same(GEN T, GEN x, GEN y)
    1508             : {
    1509     8565277 :   GEN z = cgetg(3,t_POLMOD), a;
    1510     8565277 :   long v = varn(T), lx = lg(x), ly = lg(y);
    1511     8565277 :   gel(z,1) = RgX_copy(T);
    1512             :   /* x * y mod T optimised */
    1513     8565277 :   if (typ(x) != t_POL || varn(x) != v || lx <= 3
    1514     7894312 :    || typ(y) != t_POL || varn(y) != v || ly <= 3)
    1515     1588179 :     a = gmul(x, y);
    1516             :   else
    1517             :   {
    1518     6977098 :     if (lg(T) == 5 && isint1(gel(T,4))) /* quadratic fields */
    1519     5855498 :       a = quad_polmod_mul(T, x, y);
    1520             :     else
    1521     1121600 :       a = RgXQ_mul(x, y, gel(z,1));
    1522             :   }
    1523     8565277 :   gel(z,2) = a; return z;
    1524             : }
    1525             : static GEN
    1526      484275 : sqr_polmod(GEN T, GEN x)
    1527             : {
    1528      484275 :   GEN a, z = cgetg(3,t_POLMOD);
    1529      484275 :   gel(z,1) = RgX_copy(T);
    1530      484275 :   if (typ(x) != t_POL || varn(x) != varn(T) || lg(x) <= 3)
    1531       93885 :     a = gsqr(x);
    1532             :   else
    1533             :   {
    1534      390390 :     pari_sp av = avma;
    1535      390390 :     a = RgXQ_sqr(x, gel(z,1));
    1536      390390 :     a = gerepileupto(av, a);
    1537             :   }
    1538      484275 :   gel(z,2) = a; return z;
    1539             : }
    1540             : /* Mod(x,X) * Mod(y,Y) */
    1541             : static GEN
    1542     2674917 : mul_polmod(GEN X, GEN Y, GEN x, GEN y)
    1543             : {
    1544     2674917 :   long T[3] = { evaltyp(t_POLMOD) | _evallg(3),0,0 };
    1545     2674917 :   long vx = varn(X), vy = varn(Y);
    1546     2674917 :   GEN z = cgetg(3,t_POLMOD);
    1547             : 
    1548     2674917 :   if (vx==vy) {
    1549             :     pari_sp av;
    1550          14 :     gel(z,1) = RgX_gcd(X,Y); av = avma;
    1551          14 :     warn_coercion(X,Y,gel(z,1));
    1552          14 :     gel(z,2) = gerepileupto(av, gmod(gmul(x, y), gel(z,1)));
    1553          14 :     return z;
    1554             :   }
    1555     2674903 :   if (varncmp(vx, vy) < 0)
    1556      410662 :   { gel(z,1) = RgX_copy(X); gel(T,1) = Y; gel(T,2) = y; y = T; }
    1557             :   else
    1558     2264241 :   { gel(z,1) = RgX_copy(Y); gel(T,1) = X; gel(T,2) = x; x = T; }
    1559     2674903 :   gel(z,2) = gmul(x, y); return z;
    1560             : }
    1561             : 
    1562             : #if 0 /* used by 3M only */
    1563             : /* set z = x+y and return 1 if x,y have the same sign
    1564             :  * set z = x-y and return 0 otherwise */
    1565             : static int
    1566             : did_add(GEN x, GEN y, GEN *z)
    1567             : {
    1568             :   long tx = typ(x), ty = typ(y);
    1569             :   if (tx == ty) switch(tx)
    1570             :   {
    1571             :     case t_INT: *z = addii(x,y); return 1;
    1572             :     case t_FRAC: *z = addsub_frac(x,y,addii); return 1;
    1573             :     case t_REAL:
    1574             :       if (signe(x) == -signe(y))
    1575             :       { *z = subrr(x,y); return 0; }
    1576             :       else
    1577             :       { *z = addrr(x,y); return 1; }
    1578             :   }
    1579             :   if (tx == t_REAL) switch(ty)
    1580             :   {
    1581             :     case t_INT:
    1582             :       if (signe(x) == -signe(y))
    1583             :       { *z = subri(x,y); return 0; }
    1584             :       else
    1585             :       { *z = addri(x,y); return 1; }
    1586             :     case t_FRAC:
    1587             :       if (signe(x) == -signe(gel(y,1)))
    1588             :       { *z = gsub(x,y); return 0; }
    1589             :       else
    1590             :       { *z = gadd(x,y); return 1; }
    1591             :   }
    1592             :   else if (ty == t_REAL) switch(tx)
    1593             :   {
    1594             :     case t_INT:
    1595             :       if (signe(x) == -signe(y))
    1596             :       { *z = subir(x,y); return 0; }
    1597             :       else
    1598             :       { *z = addir(x,y); return 1; }
    1599             :     case t_FRAC:
    1600             :       if (signe(gel(x,1)) == -signe(y))
    1601             :       { *z = gsub(x,y); return 0; }
    1602             :       else
    1603             :       { *z = gadd(x,y); return 1; }
    1604             :   }
    1605             :   *z = gadd(x,y); return 1;
    1606             : }
    1607             : #endif
    1608             : /* x * I * y, x t_COMPLEX with non-intzero real part, y non-intzero "scalar" */
    1609             : static GEN
    1610    11678818 : mulcIR(GEN x, GEN y)
    1611             : {
    1612    11678818 :   GEN z = cgetg(3,t_COMPLEX);
    1613    11678970 :   pari_sp av = avma;
    1614    11678970 :   gel(z,1) = gerepileupto(av, gneg(gmul(y,gel(x,2))));
    1615    11679027 :   gel(z,2) = gmul(y, gel(x,1));
    1616    11678946 :   return z;
    1617             : 
    1618             : }
    1619             : /* x,y COMPLEX */
    1620             : static GEN
    1621   273909089 : mulcc(GEN x, GEN y)
    1622             : {
    1623   273909089 :   GEN xr = gel(x,1), xi = gel(x,2);
    1624   273909089 :   GEN yr = gel(y,1), yi = gel(y,2);
    1625             :   GEN p1, p2, p3, p4, z;
    1626             :   pari_sp tetpil, av;
    1627             : 
    1628   273909089 :   if (isintzero(xr))
    1629             :   {
    1630    15657413 :     if (isintzero(yr)) {
    1631     7129092 :       av = avma;
    1632     7129092 :       return gerepileupto(av, gneg(gmul(xi,yi)));
    1633             :     }
    1634     8528155 :     return mulcIR(y, xi);
    1635             :   }
    1636   258240927 :   if (isintzero(yr)) return mulcIR(x, yi);
    1637             : 
    1638   255081780 :   z = cgetg(3,t_COMPLEX); av = avma;
    1639             : #if 0
    1640             :   /* 3M method avoiding catastrophic cancellation, BUT loses accuracy due to
    1641             :    * e.g. xr + xi if exponents differ */
    1642             :   if (did_add(xr, xi, &p3))
    1643             :   {
    1644             :     if (did_add(yr, yi, &p4)) {
    1645             :     /* R = xr*yr - xi*yi
    1646             :      * I = (xr+xi)(yr+yi) - xr*yr - xi*yi */
    1647             :       p1 = gmul(xr,yr);
    1648             :       p2 = gmul(xi,yi); p2 = gneg(p2);
    1649             :       p3 = gmul(p3, p4);
    1650             :       p4 = gsub(p2, p1);
    1651             :     } else {
    1652             :     /* R = (xr + xi) * (yr - yi) + (xr * yi - xi * yr)
    1653             :      * I = xr*yi + xi*yr */
    1654             :       p1 = gmul(p3,p4);
    1655             :       p3 = gmul(xr,yi);
    1656             :       p4 = gmul(xi,yr);
    1657             :       p2 = gsub(p3, p4);
    1658             :     }
    1659             :   } else {
    1660             :     if (did_add(yr, yi, &p4)) {
    1661             :      /* R = (xr - xi) * (yr + yi) + (xi * yr - xr * yi)
    1662             :       * I = xr*yi +xi*yr */
    1663             :       p1 = gmul(p3,p4);
    1664             :       p3 = gmul(xr,yi);
    1665             :       p4 = gmul(xi,yr);
    1666             :       p2 = gsub(p4, p3);
    1667             :     } else {
    1668             :     /* R = xr*yr - xi*yi
    1669             :      * I = -(xr-xi)(yr-yi) + xr*yr + xi*yi */
    1670             :       p3 = gneg( gmul(p3, p4) );
    1671             :       p1 = gmul(xr,yr);
    1672             :       p2 = gmul(xi,yi);
    1673             :       p4 = gadd(p1, p2);
    1674             : 
    1675             :       p2 = gneg(p2);
    1676             :     }
    1677             :   }
    1678             :   tetpil = avma;
    1679             :   gel(z,1) = gadd(p1,p2);
    1680             :   gel(z,2) = gadd(p3,p4);
    1681             : #else
    1682   255082180 :   if (typ(xr)==t_INT && typ(yr)==t_INT && typ(xi)==t_INT && typ(yi)==t_INT)
    1683             :   { /* 3M formula */
    1684      559286 :     p3 = addii(xr,xi);
    1685      559286 :     p4 = addii(yr,yi);
    1686      559286 :     p1 = mulii(xr,yr);
    1687      559286 :     p2 = mulii(xi,yi);
    1688      559286 :     p3 = mulii(p3,p4);
    1689      559286 :     p4 = addii(p2,p1);
    1690      559286 :     tetpil = avma;
    1691      559286 :     gel(z,1) = subii(p1,p2);
    1692      559286 :     gel(z,2) = subii(p3,p4);
    1693      559286 :     if (!signe(gel(z,2)))
    1694      113225 :       return gerepileuptoint((pari_sp)(z+3), gel(z,1));
    1695             :   }
    1696             :   else
    1697             :   { /* naive 4M formula: avoid all loss of accuracy */
    1698   254522894 :     p1 = gmul(xr,yr);
    1699   254465770 :     p2 = gmul(xi,yi);
    1700   254458175 :     p3 = gmul(xr,yi);
    1701   254462076 :     p4 = gmul(xi,yr);
    1702   254466345 :     tetpil = avma;
    1703   254466345 :     gel(z,1) = gsub(p1,p2);
    1704   254329685 :     gel(z,2) = gadd(p3,p4);
    1705   254336942 :     if (isintzero(gel(z,2)))
    1706             :     {
    1707       50575 :       cgiv(gel(z,2));
    1708       50575 :       return gerepileupto((pari_sp)(z+3), gel(z,1));
    1709             :     }
    1710             :   }
    1711             : #endif
    1712             : 
    1713   254718975 :   gerepilecoeffssp(av,tetpil, z+1,2); return z;
    1714             : }
    1715             : /* x,y PADIC */
    1716             : static GEN
    1717    17503317 : mulpp(GEN x, GEN y) {
    1718    17503317 :   GEN pd, p = padic_p(x), ux = padic_u(x), uy = padic_u(y);
    1719    17503317 :   long e = valp(x) + valp(y), dx, dy;
    1720             : 
    1721    17503317 :   if (!equalii(p, padic_p(y))) pari_err_OP("*",x,y);
    1722    17503238 :   if (!signe(ux) || !signe(uy)) return zeropadic(p, e);
    1723    16930392 :   dx = precp(x); dy = precp(y);
    1724    16930392 :   if (dx > dy) { pd = padic_pd(y); dx = dy; } else pd = padic_pd(x);
    1725    16930392 :   retmkpadic(Fp_mul(ux, uy, pd), icopy(p), icopy(pd), e, dx);
    1726             : }
    1727             : /* x,y QUAD */
    1728             : static GEN
    1729        1106 : mulqq(GEN x, GEN y) {
    1730        1106 :   GEN z = cgetg(4,t_QUAD);
    1731        1106 :   GEN p1, p2, p3, p4, P = gel(x,1), b = gel(P,3), c = gel(P,2);
    1732             :   pari_sp av, tetpil;
    1733        1106 :   if (!ZX_equal(P, gel(y,1))) pari_err_OP("*",x,y);
    1734             : 
    1735        1106 :   gel(z,1) = ZX_copy(P); av = avma;
    1736        1106 :   p2 = gmul(gel(x,2),gel(y,2));
    1737        1106 :   p3 = gmul(gel(x,3),gel(y,3));
    1738        1106 :   p1 = gmul(gneg_i(c),p3);
    1739             : 
    1740        1106 :   if (signe(b))
    1741         987 :     p4 = gadd(gmul(gel(x,2),gel(y,3)), gmul(gel(x,3),gel(y,2)));
    1742             :   else
    1743             :   {
    1744         119 :     p3 = gmul(gel(x,2),gel(y,3));
    1745         119 :     p4 = gmul(gel(x,3),gel(y,2));
    1746             :   }
    1747        1106 :   tetpil = avma;
    1748        1106 :   gel(z,2) = gadd(p2,p1);
    1749        1106 :   gel(z,3) = gadd(p4,p3);
    1750        1106 :   gerepilecoeffssp(av,tetpil,z+2,2); return z;
    1751             : }
    1752             : 
    1753             : GEN
    1754    15647261 : mulcxI(GEN x)
    1755             : {
    1756             :   GEN z;
    1757    15647261 :   switch(typ(x))
    1758             :   {
    1759     2632352 :     case t_INT: case t_REAL: case t_FRAC:
    1760     2632352 :       return mkcomplex(gen_0, x);
    1761    13015400 :     case t_COMPLEX:
    1762    13015400 :       if (isintzero(gel(x,1))) return gneg(gel(x,2));
    1763    13001958 :       z = cgetg(3,t_COMPLEX);
    1764    13002274 :       gel(z,1) = gneg(gel(x,2));
    1765    13002911 :       gel(z,2) = gel(x,1); return z;
    1766          48 :     default:
    1767          48 :       return gmul(gen_I(), x);
    1768             :   }
    1769             : }
    1770             : GEN
    1771     3078647 : mulcxmI(GEN x)
    1772             : {
    1773             :   GEN z;
    1774     3078647 :   switch(typ(x))
    1775             :   {
    1776      117123 :     case t_INT: case t_REAL: case t_FRAC:
    1777      117123 :       return mkcomplex(gen_0, gneg(x));
    1778     2758620 :     case t_COMPLEX:
    1779     2758620 :       if (isintzero(gel(x,1))) return gel(x,2);
    1780     2720430 :       z = cgetg(3,t_COMPLEX);
    1781     2720385 :       gel(z,1) = gel(x,2);
    1782     2720385 :       gel(z,2) = gneg(gel(x,1)); return z;
    1783      202904 :     default:
    1784      202904 :       return gmul(mkcomplex(gen_0, gen_m1), x);
    1785             :   }
    1786             : }
    1787             : /* x * I^k */
    1788             : GEN
    1789     5523936 : mulcxpowIs(GEN x, long k)
    1790             : {
    1791     5523936 :   switch (k & 3)
    1792             :   {
    1793     1369133 :     case 1: return mulcxI(x);
    1794     1362137 :     case 2: return gneg(x);
    1795     1363743 :     case 3: return mulcxmI(x);
    1796             :   }
    1797     1428923 :   return x;
    1798             : }
    1799             : 
    1800             : /* caller will assume l > 2 later */
    1801             : static GEN
    1802     5706235 : init_ser(long l, long v, long e)
    1803             : {
    1804     5706235 :   GEN z = cgetg(l, t_SER);
    1805     5706235 :   z[1] = evalvalser(e) | evalvarn(v) | evalsigne(1); return z;
    1806             : }
    1807             : 
    1808             : /* fill in coefficients of t_SER z from coeffs of t_POL y */
    1809             : static GEN
    1810     5693888 : fill_ser(GEN z, GEN y)
    1811             : {
    1812     5693888 :   long i, lx = lg(z), ly = lg(y); /* lx > 2 */
    1813     5693888 :   if (ly >= lx) {
    1814    25570501 :     for (i = 2; i < lx; i++) gel(z,i) = gel(y,i);
    1815             :   } else {
    1816      335415 :     for (i = 2; i < ly; i++) gel(z,i) = gel(y,i);
    1817      239573 :     for (     ; i < lx; i++) gel(z,i) = gen_0;
    1818             :   }
    1819             :   /* dangerous case (already normalized), don't use normalizeser */
    1820     5693888 :   if (ly == 3 && !signe(y)) { setsigne(z, 0); return z; }
    1821     5693790 :   return normalizeser(z);
    1822             : }
    1823             : /* assume typ(x) = t_VEC */
    1824             : static int
    1825          98 : is_ext_qfr(GEN x)
    1826          91 : { return lg(x) == 3 && typ(gel(x,1)) == t_QFB && !qfb_is_qfi(gel(x,1))
    1827         189 :                     && typ(gel(x,2)) == t_REAL; }
    1828             : 
    1829             : GEN
    1830  8686066884 : gmul(GEN x, GEN y)
    1831             : {
    1832             :   long tx, ty, lx, ly, vx, vy, i, l;
    1833             :   pari_sp av, tetpil;
    1834             :   GEN z, p1;
    1835             : 
    1836  8686066884 :   if (x == y) return gsqr(x);
    1837  7759276380 :   tx = typ(x); ty = typ(y);
    1838  7759276380 :   if (tx == ty) switch(tx)
    1839             :   {
    1840  3635920985 :     case t_INT: return mulii(x,y);
    1841  1869389147 :     case t_REAL: return mulrr(x,y);
    1842     1545026 :     case t_INTMOD: { GEN X = gel(x,1), Y = gel(y,1);
    1843     1545026 :       z = cgetg(3,t_INTMOD);
    1844     1545026 :       if (X==Y || equalii(X,Y))
    1845     1544991 :         return mul_intmod_same(z, X, gel(x,2), gel(y,2));
    1846          35 :       gel(z,1) = gcdii(X,Y);
    1847          35 :       warn_coercion(X,Y,gel(z,1));
    1848          35 :       av = avma; p1 = mulii(gel(x,2),gel(y,2));
    1849          35 :       gel(z,2) = gerepileuptoint(av, remii(p1, gel(z,1))); return z;
    1850             :     }
    1851    23261414 :     case t_FRAC:
    1852             :     {
    1853    23261414 :       GEN x1 = gel(x,1), x2 = gel(x,2);
    1854    23261414 :       GEN y1 = gel(y,1), y2 = gel(y,2);
    1855    23261414 :       z=cgetg(3,t_FRAC);
    1856    23261414 :       p1 = gcdii(x1, y2);
    1857    23261407 :       if (!equali1(p1)) { x1 = diviiexact(x1,p1); y2 = diviiexact(y2,p1); }
    1858    23261407 :       p1 = gcdii(x2, y1);
    1859    23261413 :       if (!equali1(p1)) { x2 = diviiexact(x2,p1); y1 = diviiexact(y1,p1); }
    1860    23261413 :       tetpil = avma;
    1861    23261413 :       gel(z,2) = mulii(x2,y2);
    1862    23261408 :       gel(z,1) = mulii(x1,y1);
    1863    23261404 :       fix_frac_if_int_GC(z,tetpil); return z;
    1864             :     }
    1865   264899242 :     case t_COMPLEX: return mulcc(x, y);
    1866    12511879 :     case t_PADIC: return mulpp(x, y);
    1867         882 :     case t_QUAD: return mulqq(x, y);
    1868    13693834 :     case t_FFELT: return FF_mul(x, y);
    1869    11004322 :     case t_POLMOD:
    1870    11004322 :       if (RgX_equal_var(gel(x,1), gel(y,1)))
    1871     8329405 :         return mul_polmod_same(gel(x,1), gel(x,2), gel(y,2));
    1872     2674917 :       return mul_polmod(gel(x,1), gel(y,1), gel(x,2), gel(y,2));
    1873   128353442 :     case t_POL:
    1874   128353442 :       vx = varn(x);
    1875   128353442 :       vy = varn(y);
    1876   128353442 :       if (vx != vy) {
    1877    45999793 :         if (varncmp(vx, vy) < 0) return RgX_Rg_mul(x, y);
    1878    26874702 :         else                     return RgX_Rg_mul(y, x);
    1879             :       }
    1880    82353649 :       return RgX_mul(x, y);
    1881             : 
    1882     4411812 :     case t_SER: {
    1883     4411812 :       vx = varn(x);
    1884     4411812 :       vy = varn(y);
    1885     4411812 :       if (vx != vy) {
    1886        3675 :         if (varncmp(vx, vy) < 0) return mul_ser_scal(x, y);
    1887        3675 :         return mul_ser_scal(y, x);
    1888             :       }
    1889     4408137 :       lx = minss(lg(x), lg(y));
    1890     4408137 :       if (lx == 2) return zeroser(vx, valser(x)+valser(y));
    1891     4408123 :       av = avma; z = init_ser(lx, vx, valser(x)+valser(y));
    1892     4408123 :       x = ser2pol_i(x, lx);
    1893     4408123 :       y = ser2pol_i(y, lx);
    1894     4408123 :       y = RgXn_mul(x, y, lx-2);
    1895     4408123 :       return gerepilecopy(av, fill_ser(z,y));
    1896             :     }
    1897          42 :     case t_VEC:
    1898          42 :       if (!is_ext_qfr(x) || !is_ext_qfr(y)) pari_err_TYPE2("*",x,y);
    1899             :     /* fall through, handle extended t_QFB */
    1900       47804 :     case t_QFB: return qfbcomp(x,y);
    1901     6720955 :     case t_RFRAC: return mul_rfrac(gel(x,1),gel(x,2), gel(y,1),gel(y,2));
    1902     4115747 :     case t_MAT: return RgM_mul(x, y);
    1903             : 
    1904        1631 :     case t_VECSMALL: /* multiply as permutation. cf perm_mul */
    1905        1631 :       z = cgetg_copy(x, &l);
    1906        1631 :       if (l != lg(y)) break;
    1907       17325 :       for (i=1; i<l; i++)
    1908             :       {
    1909       15694 :         long yi = y[i];
    1910       15694 :         if (yi < 1 || yi >= l) pari_err_TYPE2("*",x,y);
    1911       15694 :         z[i] = x[yi];
    1912             :       }
    1913        1631 :       return z;
    1914             : 
    1915           0 :     default:
    1916           0 :       pari_err_TYPE2("*",x,y);
    1917             :   }
    1918             :   /* tx != ty */
    1919  1785649032 :   if (is_const_t(ty) && is_const_t(tx))  {
    1920  1626778100 :     if (tx > ty) { swap(x,y); lswap(tx,ty); }
    1921  1626778100 :     switch(tx) {
    1922  1498107565 :     case t_INT:
    1923             :       switch(ty)
    1924             :       {
    1925   988822088 :         case t_REAL: return signe(x)? mulir(x,y): gen_0;
    1926     1290463 :         case t_INTMOD:
    1927     1290463 :           z = cgetg(3, t_INTMOD);
    1928     1290463 :           return mul_intmod_same(z, gel(y,1), gel(y,2), modii(x, gel(y,1)));
    1929    70925801 :         case t_FRAC:
    1930    70925801 :           if (!signe(x)) return gen_0;
    1931    53877003 :           z=cgetg(3,t_FRAC);
    1932    53877197 :           p1 = gcdii(x,gel(y,2));
    1933    53876323 :           if (equali1(p1))
    1934             :           {
    1935    28181425 :             set_avma((pari_sp)z);
    1936    28181414 :             gel(z,2) = icopy(gel(y,2));
    1937    28181488 :             gel(z,1) = mulii(gel(y,1), x);
    1938             :           }
    1939             :           else
    1940             :           {
    1941    25695212 :             x = diviiexact(x,p1); tetpil = avma;
    1942    25694733 :             gel(z,2) = diviiexact(gel(y,2), p1);
    1943    25694440 :             gel(z,1) = mulii(gel(y,1), x);
    1944    25694863 :             fix_frac_if_int_GC(z,tetpil);
    1945             :           }
    1946    53876867 :           return z;
    1947   434288131 :         case t_COMPLEX: return signe(x)? mulRc(x, y): gen_0;
    1948     4196901 :         case t_PADIC: return signe(x)? mulTp(x, y): gen_0;
    1949        1659 :         case t_QUAD: return mulRq(x,y);
    1950     1607167 :         case t_FFELT: return FF_Z_mul(y,x);
    1951             :       }
    1952             : 
    1953             :     case t_REAL:
    1954   121167836 :       switch(ty)
    1955             :       {
    1956    37513010 :         case t_FRAC: return mulrfrac(x, y);
    1957    83654784 :         case t_COMPLEX: return mulRc(x, y);
    1958          21 :         case t_QUAD: return mulqf(y, x, realprec(x));
    1959          21 :         default: pari_err_TYPE2("*",x,y);
    1960             :       }
    1961             : 
    1962        8036 :     case t_INTMOD:
    1963             :       switch(ty)
    1964             :       {
    1965        7133 :         case t_FRAC: { GEN X = gel(x,1);
    1966        7133 :           z = cgetg(3, t_INTMOD); p1 = Fp_mul(gel(y,1), gel(x,2), X);
    1967        7133 :           return div_intmod_same(z, X, p1, remii(gel(y,2), X));
    1968             :         }
    1969          49 :         case t_COMPLEX: return mulRc_direct(x,y);
    1970         259 :         case t_PADIC: { GEN X = gel(x,1);
    1971         259 :           z = cgetg(3, t_INTMOD);
    1972         259 :           return mul_intmod_same(z, X, gel(x,2), padic_to_Fp(y, X));
    1973             :         }
    1974          63 :         case t_QUAD: return mulRq(x, y);
    1975         532 :         case t_FFELT:
    1976         532 :           if (!equalii(gel(x,1),FF_p_i(y)))
    1977           0 :             pari_err_OP("*",x,y);
    1978         532 :           return FF_Z_mul(y,gel(x,2));
    1979             :       }
    1980             : 
    1981             :     case t_FRAC:
    1982             :       switch(ty)
    1983             :       {
    1984     5235320 :         case t_COMPLEX: return mulRc(x, y);
    1985     2582753 :         case t_PADIC: return signe(gel(x,1))? mulTp(x, y): gen_0;
    1986         343 :         case t_QUAD:  return mulRq(x, y);
    1987        2051 :         case t_FFELT: return FF_Z_Z_muldiv(y, gel(x,1),gel(x,2));
    1988             :       }
    1989             : 
    1990             :     case t_FFELT:
    1991          35 :       pari_err_TYPE2("*",x,y);
    1992             : 
    1993          21 :     case t_COMPLEX:
    1994             :       switch(ty)
    1995             :       {
    1996          14 :         case t_PADIC:
    1997          14 :           return Zp_nosquare_m1(padic_p(y))? mulRc(y, x): mulTp(x, y);
    1998           7 :         case t_QUAD:
    1999           7 :           lx = precision(x); if (!lx) pari_err_OP("*",x,y);
    2000           7 :           return mulqf(y, x, lx);
    2001             :       }
    2002             : 
    2003             :     case t_PADIC: /* ty == t_QUAD */
    2004          28 :       return (kro_quad(y,padic_p(x))== -1)? mulRq(x, y): mulTp(y, x);
    2005             :     }
    2006             :   }
    2007             : 
    2008   160538913 :   if (is_matvec_t(ty)) switch(tx)
    2009             :   {
    2010      122546 :     case t_VEC:
    2011             :       switch(ty) {
    2012       15540 :         case t_COL: return RgV_RgC_mul(x,y);
    2013      107006 :         case t_MAT: return RgV_RgM_mul(x,y);
    2014             :       }
    2015           0 :       break;
    2016        1687 :     case t_COL:
    2017             :       switch(ty) {
    2018        1687 :         case t_VEC: return RgC_RgV_mul(x,y);
    2019           0 :         case t_MAT: return RgC_RgM_mul(x,y);
    2020             :       }
    2021           0 :       break;
    2022      819646 :     case t_MAT:
    2023             :       switch(ty) {
    2024           0 :         case t_VEC: return RgM_RgV_mul(x,y);
    2025      819647 :         case t_COL: return RgM_RgC_mul(x,y);
    2026             :       }
    2027             :     default:
    2028     5402886 :       if (is_noncalc_t(tx)) pari_err_TYPE2( "*",x,y); /* necessary if ly = 1 */
    2029     5402897 :       return mul_self_scal(y, x);
    2030             :   }
    2031   157685464 :   if (is_matvec_t(tx))
    2032             :   {
    2033     2058429 :     if (is_noncalc_t(ty)) pari_err_TYPE2( "*",x,y); /* necessary if lx = 1 */
    2034     2058531 :     return mul_self_scal(x, y);
    2035             :   }
    2036   155627252 :   if (tx > ty) { swap(x,y); lswap(tx,ty); }
    2037             :   /* tx < ty, !ismatvec(x and y) */
    2038             : 
    2039   155627252 :   if (ty == t_POLMOD) /* is_const_t(tx) in this case */
    2040     3405710 :     return mul_polmod_scal(gel(y,1), gel(y,2), x);
    2041   152221542 :   if (is_scalar_t(tx)) {
    2042   148966079 :     if (tx == t_POLMOD) {
    2043     3110312 :       vx = varn(gel(x,1));
    2044     3110312 :       vy = gvar(y);
    2045     3110312 :       if (vx != vy) {
    2046     2875000 :         if (varncmp(vx,vy) > 0) return mul_scal(y, x, ty);
    2047       23093 :         return mul_polmod_scal(gel(x,1), gel(x,2), y);
    2048             :       }
    2049             :       /* error if ty == t_SER */
    2050      235312 :       av = avma; y = gmod(y, gel(x,1));
    2051      235305 :       return gerepileupto(av, mul_polmod_same(gel(x,1), gel(x,2), y));
    2052             :     }
    2053   145855767 :     return mul_scal(y, x, ty);
    2054             :   }
    2055             : 
    2056             :   /* x and y are not scalars, nor matvec */
    2057     3255307 :   vx = gvar(x);
    2058     3255320 :   vy = gvar(y);
    2059     3255320 :   if (vx != vy) /* x or y is treated as a scalar */
    2060     2781167 :     return (varncmp(vx, vy) < 0)? mul_scal(x, y, tx)
    2061     2781167 :                                 : mul_scal(y, x, ty);
    2062             :   /* vx = vy */
    2063     1835942 :   switch(tx)
    2064             :   {
    2065     1835907 :     case t_POL:
    2066             :       switch (ty)
    2067             :       {
    2068        6741 :         case t_SER:
    2069             :         {
    2070             :           long v;
    2071        6741 :           av = avma; v = RgX_valrem(x, &x);
    2072        6741 :           if (v == LONG_MAX) return gerepileupto(av, Rg_get_0(x));
    2073        6727 :           v += valser(y); ly = lg(y);
    2074        6727 :           if (ly == 2) { set_avma(av); return zeroser(vx, v); }
    2075        6496 :           if (degpol(x))
    2076             :           {
    2077        2030 :             z = init_ser(ly, vx, v);
    2078        2030 :             x = RgXn_mul(x, ser2pol_i(y, ly), ly-2);
    2079        2030 :             return gerepilecopy(av, fill_ser(z, x));
    2080             :           }
    2081             :           /* take advantage of x = c*t^v */
    2082        4466 :           set_avma(av); y = mul_ser_scal(y, gel(x,2));
    2083        4466 :           setvalser(y, v); return y;
    2084             :         }
    2085             : 
    2086     1829166 :         case t_RFRAC: return mul_rfrac_scal(gel(y,1),gel(y,2), x);
    2087             :       }
    2088           0 :       break;
    2089             : 
    2090          35 :     case t_SER:
    2091             :       switch (ty)
    2092             :       {
    2093          35 :         case t_RFRAC:
    2094          35 :           av = avma;
    2095          35 :           return gerepileupto(av, gdiv(gmul(gel(y,1),x), gel(y,2)));
    2096             :       }
    2097           0 :       break;
    2098             :   }
    2099           0 :   pari_err_TYPE2("*",x,y);
    2100             :   return NULL; /* LCOV_EXCL_LINE */
    2101             : }
    2102             : 
    2103             : /* return a nonnormalized result */
    2104             : GEN
    2105      112731 : sqr_ser_part(GEN x, long l1, long l2)
    2106             : {
    2107             :   long i, j, l;
    2108             :   pari_sp av;
    2109             :   GEN Z, z, p1, p2;
    2110             :   long mi;
    2111      112731 :   if (l2 < l1) return zeroser(varn(x), 2*valser(x));
    2112      112717 :   p2 = cgetg(l2+2, t_VECSMALL)+1; /* left on stack on exit */
    2113      112717 :   Z = cgetg(l2-l1+3, t_SER);
    2114      112717 :   Z[1] = evalvalser(2*valser(x)) | evalvarn(varn(x));
    2115      112717 :   z = Z + 2-l1;
    2116      112717 :   x += 2; mi = 0;
    2117      425390 :   for (i=0; i<l1; i++)
    2118             :   {
    2119      312673 :     p2[i] = !isrationalzero(gel(x,i)); if (p2[i]) mi = i;
    2120             :   }
    2121             : 
    2122      720212 :   for (i=l1; i<=l2; i++)
    2123             :   {
    2124      607495 :     p2[i] = !isrationalzero(gel(x,i)); if (p2[i]) mi = i;
    2125      607495 :     p1=gen_0; av=avma; l=((i+1)>>1) - 1;
    2126     1256134 :     for (j=i-mi; j<=minss(l,mi); j++)
    2127      648639 :       if (p2[j] && p2[i-j]) p1 = gadd(p1, gmul(gel(x,j),gel(x,i-j)));
    2128      607495 :     p1 = gshift(p1,1);
    2129      607495 :     if ((i&1) == 0 && p2[i>>1])
    2130       93585 :       p1 = gadd(p1, gsqr(gel(x,i>>1)));
    2131      607495 :     gel(z,i) = gerepileupto(av,p1);
    2132             :   }
    2133      112717 :   return Z;
    2134             : }
    2135             : 
    2136             : GEN
    2137  1168322447 : gsqr(GEN x)
    2138             : {
    2139             :   long i, lx;
    2140             :   pari_sp av, tetpil;
    2141             :   GEN z, p1, p2, p3, p4;
    2142             : 
    2143  1168322447 :   switch(typ(x))
    2144             :   {
    2145   955621661 :     case t_INT: return sqri(x);
    2146   195392933 :     case t_REAL: return sqrr(x);
    2147      142512 :     case t_INTMOD: { GEN X = gel(x,1);
    2148      142512 :       z = cgetg(3,t_INTMOD);
    2149      142512 :       gel(z,2) = gerepileuptoint((pari_sp)z, remii(sqri(gel(x,2)), X));
    2150      142511 :       gel(z,1) = icopy(X); return z;
    2151             :     }
    2152     4134449 :     case t_FRAC: return sqrfrac(x);
    2153             : 
    2154     8115905 :     case t_COMPLEX:
    2155     8115905 :       if (isintzero(gel(x,1))) {
    2156      238819 :         av = avma;
    2157      238819 :         return gerepileupto(av, gneg(gsqr(gel(x,2))));
    2158             :       }
    2159     7877075 :       z = cgetg(3,t_COMPLEX); av = avma;
    2160     7877074 :       p1 = gadd(gel(x,1),gel(x,2));
    2161     7876862 :       p2 = gsub(gel(x,1), gel(x,2));
    2162     7876864 :       p3 = gmul(gel(x,1),gel(x,2));
    2163     7877016 :       tetpil = avma;
    2164     7877016 :       gel(z,1) = gmul(p1,p2);
    2165     7877046 :       gel(z,2) = gshift(p3,1); gerepilecoeffssp(av,tetpil,z+1,2); return z;
    2166             : 
    2167        4795 :     case t_PADIC:
    2168             :     {
    2169        4795 :       GEN u = padic_u(x), p = padic_p(x), pd = padic_pd(x);
    2170        4795 :       long v2 = 2*valp(x), d = precp(x);
    2171        4795 :       if (!signe(u)) { x = gcopy(x); setvalp(x, v2); return x; }
    2172        4739 :       if (!absequaliu(p, 2))
    2173        3206 :         retmkpadic(Fp_sqr(u, pd), icopy(p), icopy(pd), v2, d);
    2174             :       /* p = 2*/
    2175        1533 :       if (d == 1) /* (1 + O(2))^2 = 1 + O(2^3) */
    2176           7 :         retmkpadic(gen_1, gen_2, utoipos(8), v2, 3);
    2177        1526 :       retmkpadic(remi2n(sqri(u), d + 1), gen_2, int2n(d + 1), v2, d + 1);
    2178             :     }
    2179          70 :     case t_QUAD: z = cgetg(4,t_QUAD);
    2180          70 :       p1 = gel(x,1);
    2181          70 :       gel(z,1) = ZX_copy(p1); av = avma;
    2182          70 :       p2 = gsqr(gel(x,2));
    2183          70 :       p3 = gsqr(gel(x,3));
    2184          70 :       p4 = gmul(gneg_i(gel(p1,2)),p3);
    2185             : 
    2186          70 :       if (gequal0(gel(p1,3)))
    2187             :       {
    2188           7 :         tetpil = avma;
    2189           7 :         gel(z,2) = gerepile(av,tetpil,gadd(p4,p2));
    2190           7 :         av = avma;
    2191           7 :         p2 = gmul(gel(x,2),gel(x,3)); tetpil = avma;
    2192           7 :         gel(z,3) = gerepile(av,tetpil,gmul2n(p2,1)); return z;
    2193             :       }
    2194             : 
    2195          63 :       p1 = gmul2n(gmul(gel(x,2),gel(x,3)), 1);
    2196          63 :       tetpil = avma;
    2197          63 :       gel(z,2) = gadd(p2,p4);
    2198          63 :       gel(z,3) = gadd(p1,p3);
    2199          63 :       gerepilecoeffssp(av,tetpil,z+2,2); return z;
    2200             : 
    2201      484275 :     case t_POLMOD:
    2202      484275 :       return sqr_polmod(gel(x,1), gel(x,2));
    2203             : 
    2204     2967085 :     case t_FFELT: return FF_sqr(x);
    2205             : 
    2206     1354736 :     case t_POL: return RgX_sqr(x);
    2207             : 
    2208       35014 :     case t_SER:
    2209       35014 :       lx = lg(x);
    2210       35014 :       if (ser_isexactzero(x)) {
    2211          21 :         GEN z = gcopy(x);
    2212          21 :         setvalser(z, 2*valser(x));
    2213          21 :         return z;
    2214             :       }
    2215       34993 :       if (lx < 40)
    2216       34734 :         return normalizeser( sqr_ser_part(x, 0, lx-3) );
    2217             :       else
    2218             :       {
    2219         259 :         pari_sp av = avma;
    2220         259 :         GEN z = init_ser(lx, varn(x), 2*valser(x));
    2221         259 :         x = ser2pol_i(x, lx);
    2222         259 :         x = RgXn_sqr(x, lx-2);
    2223         259 :         return gerepilecopy(av, fill_ser(z,x));
    2224             :       }
    2225             : 
    2226          63 :     case t_RFRAC: z = cgetg(3,t_RFRAC);
    2227          63 :       gel(z,1) = gsqr(gel(x,1));
    2228          63 :       gel(z,2) = gsqr(gel(x,2)); return z;
    2229             : 
    2230        1169 :     case t_MAT: return RgM_sqr(x);
    2231          28 :     case t_VEC: if (!is_ext_qfr(x)) pari_err_TYPE2("*",x,x);
    2232             :     /* fall through handle extended t_QFB */
    2233       81919 :     case t_QFB: return qfbsqr(x);
    2234         658 :     case t_VECSMALL:
    2235         658 :       z = cgetg_copy(x, &lx);
    2236       16289 :       for (i=1; i<lx; i++)
    2237             :       {
    2238       15631 :         long xi = x[i];
    2239       15631 :         if (xi < 1 || xi >= lx) pari_err_TYPE2("*",x,x);
    2240       15631 :         z[i] = x[xi];
    2241             :       }
    2242         658 :       return z;
    2243             :   }
    2244           7 :   pari_err_TYPE2("*",x,x);
    2245             :   return NULL; /* LCOV_EXCL_LINE */
    2246             : }
    2247             : 
    2248             : /********************************************************************/
    2249             : /**                                                                **/
    2250             : /**                           DIVISION                             **/
    2251             : /**                                                                **/
    2252             : /********************************************************************/
    2253             : static GEN
    2254       13162 : div_rfrac_scal(GEN x, GEN y)
    2255             : {
    2256       13162 :   pari_sp av = avma;
    2257       13162 :   GEN d = rfrac_denom_mul_scal(gel(x,2), y);
    2258       13162 :   return gerepileupto(av, gred_rfrac_simple(gel(x,1), d));
    2259             : }
    2260             : static GEN
    2261       37515 : div_scal_rfrac(GEN x, GEN y)
    2262             : {
    2263       37515 :   GEN y1 = gel(y,1), y2 = gel(y,2);
    2264       37515 :   if (typ(y1) == t_POL && varn(y2) == varn(y1))
    2265             :   {
    2266          14 :     if (degpol(y1))
    2267             :     {
    2268          14 :       pari_sp av = avma;
    2269          14 :       GEN _1 = Rg_get_1(x);
    2270          14 :       if (transtype(_1)) y1 = gmul(y1, _1);
    2271          14 :       return gerepileupto(av, gred_rfrac_simple(gmul(x, y2), y1));
    2272             :     }
    2273           0 :     y1 = gel(y1,2);
    2274             :   }
    2275       37501 :   return RgX_Rg_mul(y2, gdiv(x,y1));
    2276             : }
    2277             : static GEN
    2278     1186733 : div_rfrac(GEN x, GEN y)
    2279     1186733 : { return mul_rfrac(gel(x,1),gel(x,2), gel(y,2),gel(y,1)); }
    2280             : 
    2281             : /* x != 0 */
    2282             : static GEN
    2283     1338593 : div_ser_scal(GEN x, GEN t)
    2284             : {
    2285     1338593 :   if (ser_isexactzero(x))
    2286             :   {
    2287          28 :     GEN z = scalarser(lg(x) == 2? Rg_get_0(t): gdiv(gel(x,2), t), varn(x), 1);
    2288          28 :     setvalser(z, valser(x)); return z;
    2289             :   }
    2290     6202113 :   pari_APPLY_ser(gdiv(gel(x,i), t));
    2291             : }
    2292             : GEN
    2293        7658 : ser_normalize(GEN x)
    2294             : {
    2295        7658 :   long i, lx = lg(x);
    2296             :   GEN c, z;
    2297        7658 :   if (lx == 2) return x;
    2298        7658 :   c = gel(x,2); if (gequal1(c)) return x;
    2299        7581 :   z = cgetg(lx, t_SER); z[1] = x[1]; gel(z,2) = gen_1;
    2300      108752 :   for (i=3; i<lx; i++) gel(z,i) = gdiv(gel(x,i),c);
    2301        7581 :   return z;
    2302             : }
    2303             : 
    2304             : /* y != 0 */
    2305             : static GEN
    2306     6699017 : div_T_scal(GEN x, GEN y, long tx) {
    2307     6699017 :   switch(tx)
    2308             :   {
    2309     5349341 :     case t_POL: return RgX_Rg_div(x, y);
    2310     1338586 :     case t_SER: return div_ser_scal(x, y);
    2311       11090 :     case t_RFRAC: return div_rfrac_scal(x,y);
    2312             :   }
    2313           0 :   pari_err_TYPE2("/",x,y);
    2314             :   return NULL; /* LCOV_EXCL_LINE */
    2315             : }
    2316             : 
    2317             : static GEN
    2318     9263796 : div_scal_pol(GEN x, GEN y) {
    2319     9263796 :   long ly = lg(y);
    2320             :   pari_sp av;
    2321             :   GEN _1;
    2322     9263796 :   if (ly == 3) return scalarpol(gdiv(x,gel(y,2)), varn(y));
    2323     9184716 :   if (isrationalzero(x)) return zeropol(varn(y));
    2324     7130620 :   av = avma;
    2325     7130620 :   _1 = Rg_get_1(x); if (transtype(_1)) y = gmul(y, _1);
    2326     7130620 :   return gerepileupto(av, gred_rfrac_simple(x,y));
    2327             : }
    2328             : static GEN
    2329       18550 : div_scal_ser(GEN x, GEN y)
    2330             : {
    2331       18550 :   pari_sp av = avma;
    2332       18550 :   GEN _1 = Rg_get_1(x);
    2333       18550 :   if (transtype(_1)) y = gmul(y, _1);
    2334       18550 :   return gerepileupto(av, gmul(x, ser_inv(y)));
    2335             : }
    2336             : static GEN
    2337     9271962 : div_scal_T(GEN x, GEN y, long ty) {
    2338     9271962 :   switch(ty)
    2339             :   {
    2340     9216660 :     case t_POL: return div_scal_pol(x, y);
    2341       18543 :     case t_SER: return div_scal_ser(x, y);
    2342       36759 :     case t_RFRAC: return div_scal_rfrac(x, y);
    2343             :   }
    2344           0 :   pari_err_TYPE2("/",x,y);
    2345             :   return NULL; /* LCOV_EXCL_LINE */
    2346             : }
    2347             : 
    2348             : /* assume tx = ty = t_SER, same variable vx */
    2349             : static GEN
    2350      771510 : div_ser(GEN x, GEN y, long vx)
    2351             : {
    2352      771510 :   long e, v = valser(x) - valser(y), lx = lg(x), ly = lg(y);
    2353      771510 :   GEN y0 = y, z;
    2354      771510 :   pari_sp av = avma;
    2355             : 
    2356      771510 :   if (!signe(y)) pari_err_INV("div_ser", y);
    2357      771503 :   if (ser_isexactzero(x))
    2358             :   {
    2359       59892 :     if (lx == 2) return zeroser(vx, v);
    2360         147 :     z = scalarser(gmul(gel(x,2),Rg_get_0(y)), varn(x), 1);
    2361         147 :     setvalser(z, v); return z;
    2362             :   }
    2363      711611 :   if (lx < ly) ly = lx;
    2364      711611 :   y = ser2pol_i_normalize(y, ly, &e);
    2365      711611 :   if (e)
    2366             :   {
    2367           7 :     pari_warn(warner,"normalizing a series with 0 leading term");
    2368           7 :     ly -= e; v -= e; if (ly <= 2) pari_err_INV("div_ser", y0);
    2369             :   }
    2370      711611 :   z = init_ser(ly, vx, v);
    2371      711611 :   if (ly == 3)
    2372             :   {
    2373       12347 :     gel(z,2) = gdiv(gel(x,2), gel(y,2));
    2374       12347 :     if (gequal0(gel(z,2))) setsigne(z, 0); /* can happen: Mod(1,2)/Mod(1,3) */
    2375       12347 :     return gerepileupto(av, z);
    2376             :   }
    2377      699264 :   x = ser2pol_i(x, ly);
    2378      699264 :   y = RgXn_div_i(x, y, ly-2);
    2379      699264 :   return gerepilecopy(av, fill_ser(z,y));
    2380             : }
    2381             : /* x,y compatible PADIC */
    2382             : static GEN
    2383    13212360 : divpp(GEN x, GEN y)
    2384             : {
    2385    13212360 :   GEN M, ux = padic_u(x), uy = padic_u(y), p = padic_p(x);
    2386    13212360 :   long a = precp(x), b = precp(y), v = valp(x) - valp(y);
    2387             : 
    2388    13212360 :   if (!signe(uy)) pari_err_INV("divpp",y);
    2389    13212355 :   if (!signe(ux)) return zeropadic(p, v);
    2390    13212026 :   if (a > b) M = padic_pd(y); else { M = padic_pd(x); b = a; }
    2391    13212026 :   retmkpadic(Fp_div(ux, uy, M), icopy(padic_p(x)), icopy(M), v, b);
    2392             : }
    2393             : static GEN
    2394       37093 : div_polmod_same(GEN T, GEN x, GEN y)
    2395             : {
    2396       37093 :   long v = varn(T);
    2397       37093 :   GEN a, z = cgetg(3, t_POLMOD);
    2398       37093 :   gel(z,1) = RgX_copy(T);
    2399       37093 :   if (typ(y) != t_POL || varn(y) != v || lg(y) <= 3)
    2400       11417 :     a = gdiv(x, y);
    2401       25676 :   else if (typ(x) != t_POL || varn(x) != v || lg(x) <= 3)
    2402         119 :   {
    2403         119 :     pari_sp av = avma;
    2404         119 :     a = gerepileupto(av, gmul(x, RgXQ_inv(y, T)));
    2405             :   }
    2406       25557 :   else if (degpol(T) == 2 && isint1(gel(T,4))) /* quadratic fields */
    2407        5278 :   {
    2408        5278 :     pari_sp av = avma;
    2409        5278 :     a = quad_polmod_mul(T, x, quad_polmod_conj(y, T));
    2410        5278 :     a = RgX_Rg_div(a, quad_polmod_norm(y, T));
    2411        5278 :     a = gerepileupto(av, a);
    2412             :   }
    2413             :   else
    2414             :   {
    2415       20279 :     pari_sp av = avma;
    2416       20279 :     a = RgXQ_mul(x, ginvmod(y, gel(z,1)), gel(z,1));
    2417       20279 :     a = gerepileupto(av, a);
    2418             :   }
    2419       37093 :   gel(z,2) = a; return z;
    2420             : }
    2421             : GEN
    2422   432741037 : gdiv(GEN x, GEN y)
    2423             : {
    2424   432741037 :   long tx = typ(x), ty = typ(y), lx, ly, vx, vy, i;
    2425             :   pari_sp av, tetpil;
    2426             :   GEN z, p1, p2;
    2427             : 
    2428   432741037 :   if (tx == ty) switch(tx)
    2429             :   {
    2430    88733998 :     case t_INT:
    2431    88733998 :       return Qdivii(x,y);
    2432             : 
    2433    93661155 :     case t_REAL: return divrr(x,y);
    2434       25746 :     case t_INTMOD: { GEN X = gel(x,1), Y = gel(y,1);
    2435       25746 :       z = cgetg(3,t_INTMOD);
    2436       25746 :       if (X==Y || equalii(X,Y))
    2437       25732 :         return div_intmod_same(z, X, gel(x,2), gel(y,2));
    2438          14 :       gel(z,1) = gcdii(X,Y);
    2439          14 :       warn_coercion(X,Y,gel(z,1));
    2440          14 :       av = avma; p1 = mulii(gel(x,2), Fp_inv(gel(y,2), gel(z,1)));
    2441          14 :       gel(z,2) = gerepileuptoint(av, remii(p1, gel(z,1))); return z;
    2442             :     }
    2443     4043936 :     case t_FRAC: {
    2444     4043936 :       GEN x1 = gel(x,1), x2 = gel(x,2);
    2445     4043936 :       GEN y1 = gel(y,1), y2 = gel(y,2);
    2446     4043936 :       z = cgetg(3, t_FRAC);
    2447     4043935 :       p1 = gcdii(x1, y1);
    2448     4043933 :       if (!equali1(p1)) { x1 = diviiexact(x1,p1); y1 = diviiexact(y1,p1); }
    2449     4043933 :       p1 = gcdii(x2, y2);
    2450     4043932 :       if (!equali1(p1)) { x2 = diviiexact(x2,p1); y2 = diviiexact(y2,p1); }
    2451     4043929 :       tetpil = avma;
    2452     4043929 :       gel(z,2) = mulii(x2,y1);
    2453     4043929 :       gel(z,1) = mulii(x1,y2);
    2454     4043933 :       normalize_frac(z);
    2455     4043933 :       fix_frac_if_int_GC(z,tetpil);
    2456     4043935 :       return z;
    2457             :     }
    2458     9202053 :     case t_COMPLEX:
    2459     9202053 :       if (isintzero(gel(y,1)))
    2460             :       {
    2461      195412 :         y = gel(y,2);
    2462      195412 :         if (isintzero(gel(x,1))) return gdiv(gel(x,2), y);
    2463      116011 :         z = cgetg(3,t_COMPLEX);
    2464      116011 :         gel(z,1) = gdiv(gel(x,2), y);
    2465      116011 :         av = avma;
    2466      116011 :         gel(z,2) = gerepileupto(av, gneg(gdiv(gel(x,1), y)));
    2467      116011 :         return z;
    2468             :       }
    2469     9006635 :       av = avma; p1 = cxnorm(y); p2 = mulcc(x, conj_i(y)); tetpil = avma;
    2470     9006651 :       return gerepile(av, tetpil, gdiv(p2,p1));
    2471             : 
    2472      587462 :     case t_PADIC:
    2473      587462 :       if (!equalii(padic_p(x), padic_p(y))) pari_err_OP("/",x,y);
    2474      587462 :       return divpp(x, y);
    2475             : 
    2476         238 :     case t_QUAD:
    2477         238 :       if (!ZX_equal(gel(x,1),gel(y,1))) pari_err_OP("/",x,y);
    2478         224 :       av = avma; p1 = quadnorm(y); p2 = mulqq(x, conj_i(y)); tetpil = avma;
    2479         224 :       return gerepile(av, tetpil, gdiv(p2,p1));
    2480             : 
    2481      245916 :     case t_FFELT: return FF_div(x,y);
    2482             : 
    2483       41468 :     case t_POLMOD:
    2484       41468 :       if (RgX_equal_var(gel(x,1), gel(y,1)))
    2485       37093 :         z = div_polmod_same(gel(x,1), gel(x,2), gel(y,2));
    2486             :       else {
    2487        4375 :         av = avma;
    2488        4375 :         z = gerepileupto(av, gmul(x, ginv(y)));
    2489             :       }
    2490       41468 :       return z;
    2491             : 
    2492    18549487 :     case t_POL:
    2493    18549487 :       vx = varn(x);
    2494    18549487 :       vy = varn(y);
    2495    18549487 :       if (vx != vy) {
    2496      154994 :         if (varncmp(vx, vy) < 0) return RgX_Rg_div(x, y);
    2497       47136 :                             else return div_scal_pol(x, y);
    2498             :       }
    2499    18394493 :       if (!signe(y)) pari_err_INV("gdiv",y);
    2500    18394493 :       if (lg(y) == 3) return RgX_Rg_div(x,gel(y,2));
    2501    18266538 :       av = avma;
    2502    18266538 :       return gerepileupto(av, gred_rfrac2(x,y));
    2503             : 
    2504      771531 :     case t_SER:
    2505      771531 :       vx = varn(x);
    2506      771531 :       vy = varn(y);
    2507      771531 :       if (vx != vy) {
    2508          21 :         if (varncmp(vx, vy) < 0)
    2509             :         {
    2510          14 :           if (!signe(y)) pari_err_INV("gdiv",y);
    2511           7 :           return div_ser_scal(x, y);
    2512             :         }
    2513           7 :         return div_scal_ser(x, y);
    2514             :       }
    2515      771510 :       return div_ser(x, y, vx);
    2516     1189421 :     case t_RFRAC:
    2517     1189421 :       vx = varn(gel(x,2));
    2518     1189421 :       vy = varn(gel(y,2));
    2519     1189421 :       if (vx != vy) {
    2520        2688 :         if (varncmp(vx, vy) < 0) return div_rfrac_scal(x, y);
    2521         756 :                             else return div_scal_rfrac(x, y);
    2522             :       }
    2523     1186733 :       return div_rfrac(x,y);
    2524             : 
    2525          21 :     case t_VEC: /* handle extended t_QFB */
    2526          21 :     case t_QFB: av = avma; return gerepileupto(av, qfbcomp(x, ginv(y)));
    2527             : 
    2528          28 :     case t_MAT:
    2529          28 :       p1 = RgM_div(x,y);
    2530          28 :       if (!p1) pari_err_INV("gdiv",y);
    2531          21 :       return p1;
    2532             : 
    2533           0 :     default: pari_err_TYPE2("/",x,y);
    2534             :   }
    2535             : 
    2536   215693581 :   if (tx==t_INT && is_const_t(ty)) /* optimized for speed */
    2537             :   {
    2538     3765030 :     long s = signe(x);
    2539     3765030 :     if (!s) {
    2540      514904 :       if (gequal0(y)) pari_err_INV("gdiv",y);
    2541      514904 :       switch (ty)
    2542             :       {
    2543      511740 :       default: return gen_0;
    2544          28 :       case t_INTMOD:
    2545          28 :         z = cgetg(3,t_INTMOD);
    2546          28 :         gel(z,1) = icopy(gel(y,1));
    2547          28 :         gel(z,2) = gen_0; return z;
    2548        3136 :       case t_FFELT: return FF_zero(y);
    2549             :       }
    2550             :     }
    2551     3250126 :     if (is_pm1(x)) {
    2552     1302474 :       if (s > 0) return ginv(y);
    2553      221581 :       av = avma; return gerepileupto(av, ginv(gneg(y)));
    2554             :     }
    2555     1947652 :     switch(ty)
    2556             :     {
    2557      639385 :       case t_REAL: return divir(x,y);
    2558          21 :       case t_INTMOD:
    2559          21 :         z = cgetg(3, t_INTMOD);
    2560          21 :         return div_intmod_same(z, gel(y,1), modii(x, gel(y,1)), gel(y,2));
    2561      792256 :       case t_FRAC:
    2562      792256 :         z = cgetg(3,t_FRAC); p1 = gcdii(x,gel(y,1));
    2563      792256 :         if (equali1(p1))
    2564             :         {
    2565      271435 :           set_avma((pari_sp)z);
    2566      271435 :           gel(z,2) = icopy(gel(y,1));
    2567      271435 :           gel(z,1) = mulii(gel(y,2), x);
    2568      271435 :           normalize_frac(z);
    2569      271435 :           fix_frac_if_int(z);
    2570             :         }
    2571             :         else
    2572             :         {
    2573      520821 :           x = diviiexact(x,p1); tetpil = avma;
    2574      520821 :           gel(z,2) = diviiexact(gel(y,1), p1);
    2575      520821 :           gel(z,1) = mulii(gel(y,2), x);
    2576      520821 :           normalize_frac(z);
    2577      520821 :           fix_frac_if_int_GC(z,tetpil);
    2578             :         }
    2579      792256 :         return z;
    2580             : 
    2581         462 :       case t_FFELT: return Z_FF_div(x,y);
    2582      513797 :       case t_COMPLEX: return divRc(x,y);
    2583        1505 :       case t_PADIC: return divTp(x, y);
    2584         231 :       case t_QUAD:
    2585         231 :         av = avma; p1 = quadnorm(y); p2 = mulRq(x, conj_i(y)); tetpil = avma;
    2586         231 :         return gerepile(av, tetpil, gdiv(p2,p1));
    2587             :     }
    2588             :   }
    2589   211928545 :   if (gequal0(y))
    2590             :   {
    2591          49 :     if (is_matvec_t(tx) && lg(x) == 1) return gcopy(x);
    2592          28 :     if (ty != t_MAT) pari_err_INV("gdiv",y);
    2593             :   }
    2594             : 
    2595   211938087 :   if (is_const_t(tx) && is_const_t(ty)) switch(tx)
    2596             :   {
    2597   172606953 :     case t_REAL:
    2598   172606953 :       switch(ty)
    2599             :       {
    2600   170239533 :         case t_INT: return divri(x,y);
    2601      522503 :         case t_FRAC:
    2602      522503 :           av = avma; z = divri(mulri(x,gel(y,2)), gel(y,1));
    2603      522503 :           return gerepileuptoleaf(av, z);
    2604     1844955 :         case t_COMPLEX: return divRc(x, y);
    2605          42 :         case t_QUAD: return divfq(x, y, realprec(x));
    2606          18 :         default: pari_err_TYPE2("/",x,y);
    2607             :       }
    2608             : 
    2609         280 :     case t_INTMOD:
    2610             :       switch(ty)
    2611             :       {
    2612         196 :         case t_INT:
    2613         196 :           z = cgetg(3, t_INTMOD);
    2614         196 :           return div_intmod_same(z, gel(x,1), gel(x,2), modii(y, gel(x,1)));
    2615          28 :         case t_FRAC: { GEN X = gel(x,1);
    2616          28 :           z = cgetg(3,t_INTMOD); p1 = remii(mulii(gel(y,2), gel(x,2)), X);
    2617          28 :           return div_intmod_same(z, X, p1, modii(gel(y,1), X));
    2618             :         }
    2619           7 :         case t_FFELT:
    2620           7 :           if (!equalii(gel(x,1),FF_p_i(y)))
    2621           0 :             pari_err_OP("/",x,y);
    2622           7 :           return Z_FF_div(gel(x,2),y);
    2623             : 
    2624           0 :         case t_COMPLEX:
    2625           0 :           av = avma;
    2626           0 :           return gerepileupto(av, mulRc_direct(gdiv(x,cxnorm(y)), conj_i(y)));
    2627             : 
    2628          14 :         case t_QUAD:
    2629          14 :           av = avma; p1 = quadnorm(y); p2 = gmul(x,conj_i(y)); tetpil = avma;
    2630          14 :           return gerepile(av,tetpil, gdiv(p2,p1));
    2631             : 
    2632           7 :         case t_PADIC: { GEN X = gel(x,1);
    2633           7 :           z = cgetg(3, t_INTMOD);
    2634           7 :           return div_intmod_same(z, X, gel(x,2), padic_to_Fp(y, X));
    2635             :         }
    2636          28 :         case t_REAL: pari_err_TYPE2("/",x,y);
    2637             :       }
    2638             : 
    2639             :     case t_FRAC:
    2640             :       switch(ty)
    2641             :       {
    2642     2098255 :         case t_INT: z = cgetg(3, t_FRAC);
    2643     2098255 :         p1 = gcdii(y,gel(x,1));
    2644     2098255 :         if (equali1(p1))
    2645             :         {
    2646      828738 :           set_avma((pari_sp)z); tetpil = 0;
    2647      828738 :           gel(z,1) = icopy(gel(x,1));
    2648             :         }
    2649             :         else
    2650             :         {
    2651     1269517 :           y = diviiexact(y,p1); tetpil = avma;
    2652     1269517 :           gel(z,1) = diviiexact(gel(x,1), p1);
    2653             :         }
    2654     2098255 :         gel(z,2) = mulii(gel(x,2),y);
    2655     2098255 :         normalize_frac(z);
    2656     2098255 :         if (tetpil) fix_frac_if_int_GC(z,tetpil);
    2657     2098255 :         return z;
    2658             : 
    2659       59516 :         case t_REAL:
    2660       59516 :           av=avma; p1=mulri(y,gel(x,2)); tetpil=avma;
    2661       59516 :           return gerepile(av, tetpil, divir(gel(x,1), p1));
    2662             : 
    2663           7 :         case t_INTMOD: { GEN Y = gel(y,1);
    2664           7 :           z = cgetg(3,t_INTMOD); p1 = remii(mulii(gel(y,2),gel(x,2)), Y);
    2665           7 :           return div_intmod_same(z, Y, modii(gel(x,1), Y), p1);
    2666             :         }
    2667             : 
    2668          28 :         case t_FFELT: av=avma;
    2669          28 :           return gerepileupto(av,Z_FF_div(gel(x,1),FF_Z_mul(y,gel(x,2))));
    2670             : 
    2671         868 :         case t_COMPLEX: return divRc(x, y);
    2672             : 
    2673        2141 :         case t_PADIC:
    2674        2141 :           if (!signe(gel(x,1))) return gen_0;
    2675        2141 :           return divTp(x, y);
    2676             : 
    2677          42 :         case t_QUAD:
    2678          42 :           av=avma; p1=quadnorm(y); p2=gmul(x,conj_i(y)); tetpil=avma;
    2679          42 :           return gerepile(av,tetpil,gdiv(p2,p1));
    2680             :       }
    2681             : 
    2682             :     case t_FFELT:
    2683         161 :       switch (ty)
    2684             :       {
    2685          77 :         case t_INT: return FF_Z_Z_muldiv(x,gen_1,y);
    2686          28 :         case t_FRAC: return FF_Z_Z_muldiv(x,gel(y,2),gel(y,1));
    2687           7 :         case t_INTMOD:
    2688           7 :           if (!equalii(gel(y,1),FF_p_i(x)))
    2689           0 :             pari_err_OP("/",x,y);
    2690           7 :           return FF_Z_Z_muldiv(x,gen_1,gel(y,2));
    2691          49 :         default:
    2692          49 :         pari_err_TYPE2("/",x,y);
    2693             :       }
    2694           0 :       break;
    2695             : 
    2696    13514628 :     case t_COMPLEX:
    2697             :       switch(ty)
    2698             :       {
    2699    13514627 :         case t_INT: case t_REAL: case t_FRAC: return divcR(x,y);
    2700           0 :         case t_INTMOD: return mulRc_direct(ginv(y), x);
    2701           0 :         case t_PADIC:
    2702           0 :           return Zp_nosquare_m1(padic_p(y))? divcR(x,y): divTp(x, y);
    2703           0 :         case t_QUAD:
    2704           0 :           lx = precision(x); if (!lx) pari_err_OP("/",x,y);
    2705           0 :           return divfq(x, y, lx);
    2706             :       }
    2707             : 
    2708             :     case t_PADIC:
    2709             :       switch(ty)
    2710             :       {
    2711     1205001 :         case t_INT: case t_FRAC: { GEN p = padic_p(x);
    2712     1204910 :           return signe(padic_u(x))? divpT(x, y)
    2713     2409911 :                                   : zeropadic(p, valp(x) - Q_pval(y,p));
    2714             :         }
    2715           7 :         case t_INTMOD: { GEN Y = gel(y,1);
    2716           7 :           z = cgetg(3, t_INTMOD);
    2717           7 :           return div_intmod_same(z, Y, padic_to_Fp(x, Y), gel(y,2));
    2718             :         }
    2719          14 :         case t_COMPLEX: case t_QUAD:
    2720          14 :           av=avma; p1=gmul(x,conj_i(y)); p2=gnorm(y); tetpil=avma;
    2721          14 :           return gerepile(av,tetpil,gdiv(p1,p2));
    2722             : 
    2723          28 :         case t_REAL: pari_err_TYPE2("/",x,y);
    2724             :       }
    2725             : 
    2726             :     case t_QUAD:
    2727             :       switch (ty)
    2728             :       {
    2729        1190 :         case t_INT: case t_INTMOD: case t_FRAC:
    2730        1190 :           z = cgetg(4,t_QUAD);
    2731        1190 :           gel(z,1) = ZX_copy(gel(x,1));
    2732        1190 :           gel(z,2) = gdiv(gel(x,2), y);
    2733        1190 :           gel(z,3) = gdiv(gel(x,3), y); return z;
    2734          63 :         case t_REAL: return divqf(x, y, realprec(y));
    2735          14 :         case t_PADIC: return divTp(x, y);
    2736           0 :         case t_COMPLEX:
    2737           0 :           ly = precision(y); if (!ly) pari_err_OP("/",x,y);
    2738           0 :           return divqf(x, y, ly);
    2739             :       }
    2740             :   }
    2741    22448702 :   switch(ty) {
    2742      582344 :     case t_REAL: case t_INTMOD: case t_PADIC: case t_POLMOD:
    2743      582344 :       return gmul(x, ginv(y)); /* missing gerepile, for speed */
    2744          42 :     case t_MAT:
    2745          42 :       av = avma; p1 = RgM_inv(y);
    2746          35 :       if (!p1) pari_err_INV("gdiv",y);
    2747          28 :       return gerepileupto(av, gmul(x, p1));
    2748           0 :     case t_VEC: case t_COL:
    2749             :     case t_LIST: case t_STR: case t_VECSMALL: case t_CLOSURE:
    2750           0 :       pari_err_TYPE2("/",x,y);
    2751             :   }
    2752    21866316 :   switch(tx) {
    2753     3793859 :     case t_VEC: case t_COL: case t_MAT:
    2754    12825206 :       pari_APPLY_same(gdiv(gel(x,i),y));
    2755           0 :     case t_LIST: case t_STR: case t_VECSMALL: case t_CLOSURE:
    2756           0 :       pari_err_TYPE2("/",x,y);
    2757             :   }
    2758             : 
    2759    18072457 :   vy = gvar(y);
    2760    18073441 :   if (tx == t_POLMOD) { GEN X = gel(x,1);
    2761      209828 :     vx = varn(X);
    2762      209828 :     if (vx != vy) {
    2763      209254 :       if (varncmp(vx, vy) > 0) return div_scal_T(x, y, ty);
    2764      208841 :       retmkpolmod(gdiv(gel(x,2), y), RgX_copy(X));
    2765             :     }
    2766             :     /* y is POL, SER or RFRAC */
    2767         574 :     av = avma;
    2768         574 :     switch(ty)
    2769             :     {
    2770           0 :       case t_RFRAC: y = gmod(ginv(y), X); break;
    2771         574 :       default: y = ginvmod(gmod(y,X), X);
    2772             :     }
    2773         567 :     return gerepileupto(av, mul_polmod_same(X, gel(x,2), y));
    2774             :   }
    2775             :   /* x and y are not both is_scalar_t. If one of them is scalar, it's not a
    2776             :    * POLMOD (done already), hence its variable is NO_VARIABLE. If the other has
    2777             :    * variable NO_VARIABLE, then the operation is incorrect */
    2778    17863613 :   vx = gvar(x);
    2779    17863605 :   if (vx != vy) { /* includes cases where one is scalar */
    2780    15970577 :     if (varncmp(vx, vy) < 0) return div_T_scal(x, y, tx);
    2781     9271549 :                         else return div_scal_T(x, y, ty);
    2782             :   }
    2783     1893028 :   switch(tx)
    2784             :   {
    2785     1046171 :     case t_POL:
    2786             :       switch(ty)
    2787             :       {
    2788          28 :         case t_SER:
    2789             :         {
    2790          28 :           GEN y0 = y;
    2791             :           long v;
    2792          28 :           av = avma; v = RgX_valrem(x, &x);
    2793          28 :           if (v == LONG_MAX) return gerepileupto(av, Rg_get_0(x));
    2794          14 :           v -= valser(y); ly = lg(y); /* > 2 */
    2795          14 :           y = ser2pol_i_normalize(y, ly, &i);
    2796          14 :           if (i)
    2797             :           {
    2798           7 :             pari_warn(warner,"normalizing a series with 0 leading term");
    2799           7 :             ly -= i; v -= i; if (ly <= 2) pari_err_INV("gdiv", y0);
    2800             :           }
    2801           7 :           z = init_ser(ly, vx, v);
    2802           7 :           return gerepilecopy(av, fill_ser(z, RgXn_div(x, y, ly-2)));
    2803             :         }
    2804             : 
    2805     1046143 :         case t_RFRAC:
    2806             :         {
    2807     1046143 :           GEN y1 = gel(y,1), y2 = gel(y,2);
    2808     1046143 :           if (typ(y1) == t_POL && varn(y1) == vx)
    2809        1143 :             return mul_rfrac_scal(y2, y1, x);
    2810     1045000 :           av = avma;
    2811     1045000 :           return gerepileupto(av, RgX_Rg_div(RgX_mul(y2, x), y1));
    2812             :         }
    2813             :       }
    2814           0 :       break;
    2815             : 
    2816      584324 :     case t_SER:
    2817             :       switch(ty)
    2818             :       {
    2819      584310 :         case t_POL:
    2820             :         {
    2821      584310 :           long v = valser(x);
    2822      584310 :           lx = lg(x);
    2823      584310 :           if (lx == 2) return zeroser(vx, v - RgX_val(y));
    2824      584205 :           av = avma;
    2825      584205 :           x = ser2pol_i(x, lx); v -= RgX_valrem_inexact(y, &y);
    2826      584205 :           z = init_ser(lx, vx, v);
    2827      584205 :           if (!signe(x)) setsigne(z,0);
    2828      584205 :           return gerepilecopy(av, fill_ser(z, RgXn_div(x, y, lx - 2)));
    2829             :         }
    2830          14 :         case t_RFRAC:
    2831          14 :           av = avma;
    2832          14 :           return gerepileupto(av, gdiv(gmul(x,gel(y,2)), gel(y,1)));
    2833             :       }
    2834           0 :       break;
    2835             : 
    2836      262512 :     case t_RFRAC:
    2837             :       switch(ty)
    2838             :       {
    2839      262512 :         case t_POL: return div_rfrac_pol(gel(x,1),gel(x,2), y);
    2840           0 :         case t_SER:
    2841           0 :           av = avma;
    2842           0 :           return gerepileupto(av, gdiv(gel(x,1), gmul(gel(x,2), y)));
    2843             :       }
    2844           0 :       break;
    2845             :   }
    2846          21 :   pari_err_TYPE2("/",x,y);
    2847             :   return NULL; /* LCOV_EXCL_LINE */
    2848             : }
    2849             : 
    2850             : /********************************************************************/
    2851             : /**                                                                **/
    2852             : /**                     SIMPLE MULTIPLICATION                      **/
    2853             : /**                                                                **/
    2854             : /********************************************************************/
    2855             : GEN
    2856   230008921 : gmulsg(long s, GEN x)
    2857             : {
    2858             :   pari_sp av;
    2859             :   long i;
    2860             :   GEN z;
    2861             : 
    2862   230008921 :   switch(typ(x))
    2863             :   {
    2864   137227618 :     case t_INT:  return mulsi(s,x);
    2865    69457463 :     case t_REAL: return s? mulsr(s,x): gen_0; /* gmul semantic */
    2866      174411 :     case t_INTMOD: { GEN p = gel(x,1);
    2867      174411 :       z = cgetg(3,t_INTMOD);
    2868      174411 :       gel(z,2) = gerepileuptoint((pari_sp)z, modii(mulsi(s,gel(x,2)), p));
    2869      174409 :       gel(z,1) = icopy(p); return z;
    2870             :     }
    2871      548083 :     case t_FFELT: return FF_Z_mul(x,stoi(s));
    2872     6210802 :     case t_FRAC:
    2873     6210802 :       if (!s) return gen_0;
    2874     6196102 :       z = cgetg(3,t_FRAC);
    2875     6196102 :       i = labs(s); i = ugcdiu(gel(x,2), i);
    2876     6196102 :       if (i == 1)
    2877             :       {
    2878     2269533 :         gel(z,2) = icopy(gel(x,2));
    2879     2269533 :         gel(z,1) = mulis(gel(x,1), s);
    2880             :       }
    2881             :       else
    2882             :       {
    2883     3926569 :         gel(z,2) = diviuexact(gel(x,2), (ulong)i);
    2884     3926569 :         gel(z,1) = mulis(gel(x,1), s/i);
    2885     3926569 :         fix_frac_if_int(z);
    2886             :       }
    2887     6196102 :       return z;
    2888             : 
    2889    14211681 :     case t_COMPLEX:
    2890    14211681 :       if (!s) return gen_0;
    2891    13085877 :       z = cgetg(3, t_COMPLEX);
    2892    13085859 :       gel(z,1) = gmulsg(s,gel(x,1));
    2893    13084380 :       gel(z,2) = gmulsg(s,gel(x,2)); return z;
    2894             : 
    2895        1449 :     case t_PADIC:
    2896        1449 :       if (!s) return gen_0;
    2897        1449 :       av = avma; return gerepileupto(av, mulpp(cvtop2(stoi(s),x), x));
    2898             : 
    2899           7 :     case t_QUAD: z = cgetg(4, t_QUAD);
    2900           7 :       gel(z,1) = ZX_copy(gel(x,1));
    2901           7 :       gel(z,2) = gmulsg(s,gel(x,2));
    2902           7 :       gel(z,3) = gmulsg(s,gel(x,3)); return z;
    2903             : 
    2904      204442 :     case t_POLMOD:
    2905      204442 :       retmkpolmod(gmulsg(s,gel(x,2)), RgX_copy(gel(x,1)));
    2906             : 
    2907      748902 :     case t_POL:
    2908      748902 :       if (!signe(x)) return RgX_copy(x);
    2909      732102 :       if (!s) return scalarpol(Rg_get_0(x), varn(x));
    2910     3096270 :       pari_APPLY_pol(gmulsg(s,gel(x,i)));
    2911             : 
    2912         182 :     case t_SER:
    2913         182 :       if (ser_isexactzero(x)) return gcopy(x);
    2914         182 :       if (!s) return Rg_get_0(x);
    2915        3864 :       pari_APPLY_ser(gmulsg(s,gel(x,i)));
    2916             : 
    2917           0 :     case t_RFRAC:
    2918           0 :       if (!s) return zeropol(varn(gel(x,2)));
    2919           0 :       if (s == 1) return gcopy(x);
    2920           0 :       if (s == -1) return gneg(x);
    2921           0 :       return mul_rfrac_scal(gel(x,1), gel(x,2), stoi(s));
    2922             : 
    2923     1234133 :     case t_VEC: case t_COL: case t_MAT:
    2924     3898052 :       pari_APPLY_same(gmulsg(s,gel(x,i)));
    2925             :   }
    2926           0 :   pari_err_TYPE("gmulsg",x);
    2927             :   return NULL; /* LCOV_EXCL_LINE */
    2928             : }
    2929             : 
    2930             : GEN
    2931   124340541 : gmulug(ulong s, GEN x)
    2932             : {
    2933             :   pari_sp av;
    2934             :   long i;
    2935             :   GEN z;
    2936             : 
    2937   124340541 :   switch(typ(x))
    2938             :   {
    2939   122223407 :     case t_INT:  return mului(s,x);
    2940     1059371 :     case t_REAL: return s? mulur(s,x): gen_0; /* gmul semantic */
    2941         364 :     case t_INTMOD: { GEN p = gel(x,1);
    2942         364 :       z = cgetg(3,t_INTMOD);
    2943         364 :       gel(z,2) = gerepileuptoint((pari_sp)z, modii(mului(s,gel(x,2)), p));
    2944         364 :       gel(z,1) = icopy(p); return z;
    2945             :     }
    2946         413 :     case t_FFELT: return FF_Z_mul(x,utoi(s));
    2947      981585 :     case t_FRAC:
    2948      981585 :       if (!s) return gen_0;
    2949      981571 :       z = cgetg(3,t_FRAC);
    2950      981571 :       i = ugcdiu(gel(x,2), s);
    2951      981571 :       if (i == 1)
    2952             :       {
    2953      822559 :         gel(z,2) = icopy(gel(x,2));
    2954      822559 :         gel(z,1) = muliu(gel(x,1), s);
    2955             :       }
    2956             :       else
    2957             :       {
    2958      159012 :         gel(z,2) = diviuexact(gel(x,2), i);
    2959      159012 :         gel(z,1) = muliu(gel(x,1), s/i);
    2960      159012 :         fix_frac_if_int(z);
    2961             :       }
    2962      981571 :       return z;
    2963             : 
    2964       30765 :     case t_COMPLEX:
    2965       30765 :       if (!s) return gen_0;
    2966       30765 :       z = cgetg(3, t_COMPLEX);
    2967       30765 :       gel(z,1) = gmulug(s,gel(x,1));
    2968       30765 :       gel(z,2) = gmulug(s,gel(x,2)); return z;
    2969             : 
    2970        7342 :     case t_PADIC:
    2971        7342 :       if (!s) return gen_0;
    2972        7342 :       av = avma; return gerepileupto(av, mulpp(cvtop2(utoi(s),x), x));
    2973             : 
    2974           0 :     case t_QUAD: z = cgetg(4, t_QUAD);
    2975           0 :       gel(z,1) = ZX_copy(gel(x,1));
    2976           0 :       gel(z,2) = gmulug(s,gel(x,2));
    2977           0 :       gel(z,3) = gmulug(s,gel(x,3)); return z;
    2978             : 
    2979        6783 :     case t_POLMOD:
    2980        6783 :       retmkpolmod(gmulug(s,gel(x,2)), RgX_copy(gel(x,1)));
    2981             : 
    2982       16226 :     case t_POL:
    2983       16226 :       if (!signe(x)) return RgX_copy(x);
    2984       15449 :       if (!s) return scalarpol(Rg_get_0(x), varn(x));
    2985       52843 :       pari_APPLY_pol(gmulug(s,gel(x,i)));
    2986             : 
    2987           0 :     case t_SER:
    2988           0 :       if (ser_isexactzero(x)) return gcopy(x);
    2989           0 :       if (!s) return Rg_get_0(x);
    2990           0 :       pari_APPLY_ser(gmulug(s,gel(x,i)));
    2991             : 
    2992           0 :     case t_RFRAC:
    2993           0 :       if (!s) return zeropol(varn(gel(x,2)));
    2994           0 :       if (s == 1) return gcopy(x);
    2995           0 :       return mul_rfrac_scal(gel(x,1), gel(x,2), utoi(s));
    2996             : 
    2997       14286 :     case t_VEC: case t_COL: case t_MAT:
    2998       74459 :       pari_APPLY_same(gmulug(s,gel(x,i)));
    2999             :   }
    3000           0 :   pari_err_TYPE("gmulsg",x);
    3001             :   return NULL; /* LCOV_EXCL_LINE */
    3002             : }
    3003             : 
    3004             : /********************************************************************/
    3005             : /**                                                                **/
    3006             : /**                       SIMPLE DIVISION                          **/
    3007             : /**                                                                **/
    3008             : /********************************************************************/
    3009             : 
    3010             : GEN
    3011    12939131 : gdivgs(GEN x, long s)
    3012             : {
    3013    12939131 :   long tx = typ(x), i;
    3014             :   pari_sp av;
    3015             :   GEN z;
    3016             : 
    3017    12939131 :   if (!s)
    3018             :   {
    3019           0 :     if (is_matvec_t(tx) && lg(x) == 1) return gcopy(x);
    3020           0 :     pari_err_INV("gdivgs",gen_0);
    3021             :   }
    3022    12939138 :   switch(tx)
    3023             :   {
    3024     1561457 :     case t_INT: return Qdivis(x, s);
    3025     8499281 :     case t_REAL: return divrs(x,s);
    3026             : 
    3027         357 :     case t_INTMOD:
    3028         357 :       z = cgetg(3, t_INTMOD);
    3029         357 :       return div_intmod_same(z, gel(x,1), gel(x,2), modsi(s, gel(x,1)));
    3030             : 
    3031         735 :     case t_FFELT: return FF_Z_Z_muldiv(x,gen_1,stoi(s));
    3032             : 
    3033      552150 :     case t_FRAC: z = cgetg(3, t_FRAC);
    3034      552150 :       i = ugcdiu(gel(x,1), labs(s));
    3035      552150 :       if (i == 1)
    3036             :       {
    3037      409320 :         gel(z,2) = mulsi(s, gel(x,2));
    3038      409320 :         gel(z,1) = icopy(gel(x,1));
    3039             :       }
    3040             :       else
    3041             :       {
    3042      142830 :         gel(z,2) = mulsi(s/i, gel(x,2));
    3043      142830 :         gel(z,1) = divis(gel(x,1), i);
    3044             :       }
    3045      552150 :       normalize_frac(z);
    3046      552150 :       fix_frac_if_int(z); return z;
    3047             : 
    3048     1841734 :     case t_COMPLEX: z = cgetg(3, t_COMPLEX);
    3049     1841734 :       gel(z,1) = gdivgs(gel(x,1),s);
    3050     1841734 :       gel(z,2) = gdivgs(gel(x,2),s); return z;
    3051             : 
    3052         133 :     case t_PADIC: /* divpT */
    3053             :     {
    3054         133 :       GEN p = padic_p(x);
    3055         133 :       if (!signe(padic_u(x))) return zeropadic(p, valp(x) - u_pval(s,p));
    3056         133 :       av = avma;
    3057         133 :       return gerepileupto(av, divpp(x, cvtop2(stoi(s),x)));
    3058             :     }
    3059             : 
    3060          28 :     case t_QUAD: z = cgetg(4, t_QUAD);
    3061          28 :       gel(z,1) = ZX_copy(gel(x,1));
    3062          28 :       gel(z,2) = gdivgs(gel(x,2),s);
    3063          28 :       gel(z,3) = gdivgs(gel(x,3),s); return z;
    3064             : 
    3065       37891 :     case t_POLMOD:
    3066       37891 :       retmkpolmod(gdivgs(gel(x,2),s), RgX_copy(gel(x,1)));
    3067             : 
    3068          91 :     case t_RFRAC:
    3069          91 :       if (s == 1) return gcopy(x);
    3070          84 :       else if (s == -1) return gneg(x);
    3071          84 :       return div_rfrac_scal(x, stoi(s));
    3072             : 
    3073      157686 :     case t_POL: pari_APPLY_pol_normalized(gdivgs(gel(x,i),s));
    3074           0 :     case t_SER: pari_APPLY_ser_normalized(gdivgs(gel(x,i),s));
    3075      407593 :     case t_VEC:
    3076             :     case t_COL:
    3077      894977 :     case t_MAT: pari_APPLY_same(gdivgs(gel(x,i),s));
    3078             :   }
    3079           0 :   pari_err_TYPE2("/",x, stoi(s));
    3080             :   return NULL; /* LCOV_EXCL_LINE */
    3081             : }
    3082             : 
    3083             : GEN
    3084    54715893 : gdivgu(GEN x, ulong s)
    3085             : {
    3086    54715893 :   long tx = typ(x), i;
    3087             :   pari_sp av;
    3088             :   GEN z;
    3089             : 
    3090    54715893 :   if (!s)
    3091             :   {
    3092           0 :     if (is_matvec_t(tx) && lg(x) == 1) return gcopy(x);
    3093           0 :     pari_err_INV("gdivgu",gen_0);
    3094             :   }
    3095    54715806 :   switch(tx)
    3096             :   {
    3097    17916105 :     case t_INT: return Qdiviu(x, s);
    3098    12956490 :     case t_REAL: return divru(x,s);
    3099             : 
    3100      210315 :     case t_INTMOD:
    3101      210315 :       z = cgetg(3, t_INTMOD); s = umodui(s, gel(x,1));
    3102      210315 :       return div_intmod_same(z, gel(x,1), gel(x,2), utoi(s));
    3103             : 
    3104         308 :     case t_FFELT: return FF_Z_Z_muldiv(x,gen_1,utoi(s));
    3105             : 
    3106      642010 :     case t_FRAC: z = cgetg(3, t_FRAC);
    3107      642010 :       i = ugcdiu(gel(x,1), s);
    3108      642010 :       if (i == 1)
    3109             :       {
    3110      462134 :         gel(z,2) = mului(s, gel(x,2));
    3111      462134 :         gel(z,1) = icopy(gel(x,1));
    3112             :       }
    3113             :       else
    3114             :       {
    3115      179876 :         gel(z,2) = mului(s/i, gel(x,2));
    3116      179876 :         gel(z,1) = divis(gel(x,1), i);
    3117             :       }
    3118      642010 :       normalize_frac(z);
    3119      642010 :       fix_frac_if_int(z); return z;
    3120             : 
    3121    11384462 :     case t_COMPLEX: z = cgetg(3, t_COMPLEX);
    3122    11384462 :       gel(z,1) = gdivgu(gel(x,1),s);
    3123    11384462 :       gel(z,2) = gdivgu(gel(x,2),s); return z;
    3124             : 
    3125    11551411 :     case t_PADIC: /* divpT */
    3126             :     {
    3127    11551411 :       GEN p = padic_p(x);
    3128    11551411 :       if (!signe(padic_u(x))) return zeropadic(p, valp(x) - u_pval(s,p));
    3129    11415796 :       av = avma;
    3130    11415796 :       return gerepileupto(av, divpp(x, cvtop2(utoi(s),x)));
    3131             :     }
    3132             : 
    3133           0 :     case t_QUAD: z = cgetg(4, t_QUAD);
    3134           0 :       gel(z,1) = ZX_copy(gel(x,1));
    3135           0 :       gel(z,2) = gdivgu(gel(x,2),s);
    3136           0 :       gel(z,3) = gdivgu(gel(x,3),s); return z;
    3137             : 
    3138        1456 :     case t_POLMOD:
    3139        1456 :       retmkpolmod(gdivgu(gel(x,2),s), RgX_copy(gel(x,1)));
    3140             : 
    3141          56 :     case t_RFRAC:
    3142          56 :       if (s == 1) return gcopy(x);
    3143          56 :       return div_rfrac_scal(x, utoi(s));
    3144             : 
    3145      205660 :     case t_POL: pari_APPLY_pol_normalized(gdivgu(gel(x,i),s));
    3146       16730 :     case t_SER: pari_APPLY_ser_normalized(gdivgu(gel(x,i),s));
    3147         329 :     case t_VEC:
    3148             :     case t_COL:
    3149        1148 :     case t_MAT: pari_APPLY_same(gdivgu(gel(x,i),s));
    3150             :   }
    3151           0 :   pari_err_TYPE2("/",x, utoi(s));
    3152             :   return NULL; /* LCOV_EXCL_LINE */
    3153             : }
    3154             : 
    3155             : /* x / (i*(i+1)) */
    3156             : GEN
    3157   224326065 : divrunextu(GEN x, ulong i)
    3158             : {
    3159   224326065 :   if (i & HIGHMASK) /* i(i+1) >= 2^BITS_IN_LONG*/
    3160           0 :     return divri(x, muluu(i , i+1));
    3161             :   else
    3162   224326065 :     return divru(x, i*(i+1));
    3163             : }
    3164             : /* x / (i*(i+1)) */
    3165             : GEN
    3166      814832 : gdivgunextu(GEN x, ulong i)
    3167             : {
    3168      814832 :   if (i & HIGHMASK) /* i(i+1) >= 2^BITS_IN_LONG*/
    3169           0 :     return gdivgu(x, i*(i+1));
    3170             :   else
    3171      814832 :     return gdiv(x, muluu(i, i+1));
    3172             : }
    3173             : 
    3174             : /* True shift (exact multiplication by 2^n) */
    3175             : GEN
    3176   133238500 : gmul2n(GEN x, long n)
    3177             : {
    3178             :   GEN z, a, b;
    3179             :   long k, l;
    3180             : 
    3181   133238500 :   switch(typ(x))
    3182             :   {
    3183    36365542 :     case t_INT:
    3184    36365542 :       if (n>=0) return shifti(x,n);
    3185     5051940 :       if (!signe(x)) return gen_0;
    3186     3353179 :       l = vali(x); n = -n;
    3187     3353240 :       if (n<=l) return shifti(x,-n);
    3188      398381 :       z = cgetg(3,t_FRAC);
    3189      398381 :       gel(z,1) = shifti(x,-l);
    3190      398381 :       gel(z,2) = int2n(n-l); return z;
    3191             : 
    3192    62833199 :     case t_REAL:
    3193    62833199 :       return shiftr(x,n);
    3194             : 
    3195      181023 :     case t_INTMOD: b = gel(x,1); a = gel(x,2);
    3196      181023 :       z = cgetg(3,t_INTMOD);
    3197      181023 :       if (n <= 0) return div_intmod_same(z, b, a, modii(int2n(-n), b));
    3198       82862 :       gel(z,2) = gerepileuptoint((pari_sp)z, modii(shifti(a,n), b));
    3199       82862 :       gel(z,1) = icopy(b); return z;
    3200             : 
    3201      217419 :     case t_FFELT: return FF_mul2n(x,n);
    3202             : 
    3203     1271240 :     case t_FRAC: a = gel(x,1); b = gel(x,2);
    3204     1271240 :       l = vali(a);
    3205     1271240 :       k = vali(b);
    3206     1271240 :       if (n+l >= k)
    3207             :       {
    3208      396723 :         if (expi(b) == k) return shifti(a,n-k); /* b power of 2 */
    3209      303934 :         l = n-k; k = -k;
    3210             :       }
    3211             :       else
    3212             :       {
    3213      874517 :         k = -(l+n); l = -l;
    3214             :       }
    3215     1178451 :       z = cgetg(3,t_FRAC);
    3216     1178451 :       gel(z,1) = shifti(a,l);
    3217     1178451 :       gel(z,2) = shifti(b,k); return z;
    3218             : 
    3219    10727178 :     case t_COMPLEX: z = cgetg(3,t_COMPLEX);
    3220    10727179 :       gel(z,1) = gmul2n(gel(x,1),n);
    3221    10727180 :       gel(z,2) = gmul2n(gel(x,2),n); return z;
    3222             : 
    3223         105 :     case t_QUAD: z = cgetg(4,t_QUAD);
    3224         105 :       gel(z,1) = ZX_copy(gel(x,1));
    3225         105 :       gel(z,2) = gmul2n(gel(x,2),n);
    3226         105 :       gel(z,3) = gmul2n(gel(x,3),n); return z;
    3227             : 
    3228      193830 :     case t_POLMOD:
    3229      193830 :       retmkpolmod(gmul2n(gel(x,2),n), RgX_copy(gel(x,1)));
    3230             : 
    3231     5063150 :     case t_POL:
    3232    21510330 :       pari_APPLY_pol(gmul2n(gel(x,i),n));
    3233      102847 :     case t_SER:
    3234      102847 :       if (ser_isexactzero(x)) return gcopy(x);
    3235      634256 :       pari_APPLY_ser(gmul2n(gel(x,i),n));
    3236    16279337 :     case t_VEC: case t_COL: case t_MAT:
    3237    60002133 :       pari_APPLY_same(gmul2n(gel(x,i),n));
    3238             : 
    3239          21 :     case t_RFRAC: /* int2n wrong if n < 0 */
    3240          21 :       return mul_rfrac_scal(gel(x,1),gel(x,2), gmul2n(gen_1,n));
    3241             : 
    3242        5628 :     case t_PADIC: /* int2n wrong if n < 0 */
    3243        5628 :       return gmul(gmul2n(gen_1,n),x);
    3244             :   }
    3245           0 :   pari_err_TYPE("gmul2n",x);
    3246             :   return NULL; /* LCOV_EXCL_LINE */
    3247             : }
    3248             : 
    3249             : /*******************************************************************/
    3250             : /*                                                                 */
    3251             : /*                              INVERSE                            */
    3252             : /*                                                                 */
    3253             : /*******************************************************************/
    3254             : static GEN
    3255      209970 : inv_polmod(GEN T, GEN x)
    3256             : {
    3257      209970 :   GEN z = cgetg(3,t_POLMOD), a;
    3258      209970 :   gel(z,1) = RgX_copy(T);
    3259      209969 :   if (typ(x) != t_POL || varn(x) != varn(T) || lg(x) <= 3)
    3260       85329 :     a = ginv(x);
    3261             :   else
    3262             :   {
    3263      124640 :     if (lg(T) == 5) /* quadratic fields */
    3264       13062 :       a = RgX_Rg_div(quad_polmod_conj(x,T), quad_polmod_norm(x,T));
    3265             :     else
    3266      111578 :       a = RgXQ_inv(x, T);
    3267             :   }
    3268      209969 :   gel(z,2) = a; return z;
    3269             : }
    3270             : 
    3271             : GEN
    3272    36379314 : ginv(GEN x)
    3273             : {
    3274             :   long s;
    3275             :   pari_sp av, tetpil;
    3276             :   GEN z, y, p1, p2;
    3277             : 
    3278    36379314 :   switch(typ(x))
    3279             :   {
    3280     9404859 :     case t_INT:
    3281     9404859 :       if (is_pm1(x)) return icopy(x);
    3282     7883324 :       s = signe(x); if (!s) pari_err_INV("ginv",gen_0);
    3283     7883310 :       z = cgetg(3,t_FRAC);
    3284     7889319 :       gel(z,1) = s<0? gen_m1: gen_1;
    3285     7889319 :       gel(z,2) = absi(x); return z;
    3286             : 
    3287    10751847 :     case t_REAL: return invr(x);
    3288             : 
    3289        6062 :     case t_INTMOD: z=cgetg(3,t_INTMOD);
    3290        6062 :       gel(z,1) = icopy(gel(x,1));
    3291        6062 :       gel(z,2) = Fp_inv(gel(x,2),gel(x,1)); return z;
    3292             : 
    3293      532063 :     case t_FRAC: {
    3294      532063 :       GEN a = gel(x,1), b = gel(x,2);
    3295      532063 :       s = signe(a);
    3296      532063 :       if (is_pm1(a)) return s > 0? icopy(b): negi(b);
    3297      247539 :       z = cgetg(3,t_FRAC);
    3298      247540 :       gel(z,1) = icopy(b);
    3299      247539 :       gel(z,2) = icopy(a);
    3300      247540 :       normalize_frac(z); return z;
    3301             :     }
    3302     9612947 :     case t_COMPLEX:
    3303     9612947 :       av=avma;
    3304     9612947 :       p1=cxnorm(x);
    3305     9612041 :       p2=mkcomplex(gel(x,1), gneg(gel(x,2)));
    3306     9612808 :       tetpil=avma;
    3307     9612808 :       return gerepile(av,tetpil,divcR(p2,p1));
    3308             : 
    3309         273 :     case t_QUAD:
    3310         273 :       av=avma; p1=quadnorm(x); p2=conj_i(x); tetpil=avma;
    3311         273 :       return gerepile(av,tetpil,gdiv(p2,p1));
    3312             : 
    3313      358358 :     case t_PADIC:
    3314             :     {
    3315      358358 :       GEN p = padic_p(x), pd = padic_pd(x), u = padic_u(x);
    3316      358358 :       long d = precp(x);
    3317      358358 :       if (!signe(u)) pari_err_INV("ginv",x);
    3318      358351 :       retmkpadic(Zp_inv(u, p, d), icopy(p), icopy(pd), -valp(x), d);
    3319             :     }
    3320             : 
    3321      209970 :     case t_POLMOD: return inv_polmod(gel(x,1), gel(x,2));
    3322       14027 :     case t_FFELT: return FF_inv(x);
    3323     5322884 :     case t_POL: return gred_rfrac_simple(gen_1,x);
    3324       26341 :     case t_SER: return ser_inv(x);
    3325        2933 :     case t_RFRAC:
    3326             :     {
    3327        2933 :       GEN n = gel(x,1), d = gel(x,2);
    3328        2933 :       pari_sp av = avma, ltop;
    3329        2933 :       if (gequal0(n)) pari_err_INV("ginv",x);
    3330             : 
    3331        2933 :       n = simplify_shallow(n);
    3332        2933 :       if (typ(n) != t_POL || varn(n) != varn(d))
    3333             :       {
    3334        2933 :         if (gequal1(n)) { set_avma(av); return RgX_copy(d); }
    3335         679 :         ltop = avma;
    3336         679 :         z = RgX_Rg_div(d,n);
    3337             :       } else {
    3338           0 :         ltop = avma;
    3339           0 :         z = cgetg(3,t_RFRAC);
    3340           0 :         gel(z,1) = RgX_copy(d);
    3341           0 :         gel(z,2) = RgX_copy(n);
    3342             :       }
    3343         679 :       stackdummy(av, ltop);
    3344         679 :       return z;
    3345             :     }
    3346             : 
    3347          21 :     case t_VEC: if (!is_ext_qfr(x)) break;
    3348             :     case t_QFB:
    3349       99588 :       return qfbpow(x, gen_m1);
    3350       38837 :     case t_MAT:
    3351       38837 :       y = RgM_inv(x);
    3352       38830 :       if (!y) pari_err_INV("ginv",x);
    3353       38760 :       return y;
    3354          28 :     case t_VECSMALL:
    3355             :     {
    3356          28 :       long i, lx = lg(x)-1;
    3357          28 :       y = zero_zv(lx);
    3358         112 :       for (i=1; i<=lx; i++)
    3359             :       {
    3360          84 :         long xi = x[i];
    3361          84 :         if (xi<1 || xi>lx || y[xi])
    3362           0 :           pari_err_TYPE("ginv [not a permutation]", x);
    3363          84 :         y[xi] = i;
    3364             :       }
    3365          28 :       return y;
    3366             :     }
    3367             :   }
    3368           6 :   pari_err_TYPE("inverse",x);
    3369             :   return NULL; /* LCOV_EXCL_LINE */
    3370             : }

Generated by: LCOV version 1.16