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 - gen2.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.18.0 lcov report (development 29796-c9bd12ad19) Lines: 1568 1780 88.1 %
Date: 2024-12-09 09:09:46 Functions: 148 159 93.1 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation; 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             : /**                        (second part)                           **/
      19             : /**                                                                **/
      20             : /********************************************************************/
      21             : #include "pari.h"
      22             : #include "paripriv.h"
      23             : 
      24             : /*********************************************************************/
      25             : /**                                                                 **/
      26             : /**                MAP FUNCTIONS WITH GIVEN PROTOTYPES              **/
      27             : /**                                                                 **/
      28             : /*********************************************************************/
      29             : GEN
      30         462 : map_proto_G(GEN (*f)(GEN), GEN x)
      31             : {
      32         462 :   if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_G(f, gel(x,i)));
      33         462 :   return f(x);
      34             : }
      35             : 
      36             : GEN
      37    37180571 : map_proto_lG(long (*f)(GEN), GEN x)
      38             : {
      39    37180655 :   if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_lG(f, gel(x,i)));
      40    37174137 :   return stoi(f(x));
      41             : }
      42             : 
      43             : GEN
      44         126 : map_proto_lGL(long (*f)(GEN,long), GEN x, long y)
      45             : {
      46         154 :   if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_lGL(f,gel(x,i),y));
      47         119 :   return stoi(f(x,y));
      48             : }
      49             : 
      50             : static GEN
      51     2205255 : _domul(void *data, GEN x, GEN y)
      52             : {
      53     2205255 :   GEN (*mul)(GEN,GEN)=(GEN (*)(GEN,GEN)) data;
      54     2205255 :   return mul(x,y);
      55             : }
      56             : 
      57             : GEN
      58     2419015 : gassoc_proto(GEN (*f)(GEN,GEN), GEN x, GEN y)
      59             : {
      60     2419015 :   if (!y)
      61             :   {
      62     2419015 :     pari_sp av = avma;
      63     2419015 :     switch(typ(x))
      64             :     {
      65          21 :       case t_LIST:
      66          21 :         x = list_data(x); if (!x) return gen_1;
      67             :       case t_VEC:
      68     2419001 :       case t_COL: break;
      69           7 :       default: pari_err_TYPE("association",x);
      70             :     }
      71     2419001 :     return gerepileupto(av, gen_product(x, (void *)f, _domul));
      72             : 
      73             :   }
      74           0 :   return f(x,y);
      75             : }
      76             : /*******************************************************************/
      77             : /*                                                                 */
      78             : /*                    CREATION OF A P-ADIC GEN                     */
      79             : /*                                                                 */
      80             : /*******************************************************************/
      81             : GEN
      82    16931778 : cgetp(GEN x)
      83             : {
      84    16931778 :   GEN y = cgetg(5,t_PADIC);
      85    16931753 :   y[1] = (x[1]&PRECPBITS) | _evalvalp(0);
      86    16931753 :   gel(y,2) = icopy(gel(x,2));
      87    16931697 :   gel(y,3) = icopy(gel(x,3));
      88    16931705 :   gel(y,4) = cgeti(lgefint(gel(x,3))); return y;
      89             : }
      90             : 
      91             : /*******************************************************************/
      92             : /*                                                                 */
      93             : /*                            SIZES                                */
      94             : /*                                                                 */
      95             : /*******************************************************************/
      96             : 
      97             : long
      98     5145627 : glength(GEN x)
      99             : {
     100     5145627 :   long tx = typ(x);
     101     5145627 :   switch(tx)
     102             :   {
     103         126 :     case t_INT:  return lgefint(x)-2;
     104         609 :     case t_LIST: {
     105         609 :       GEN L = list_data(x);
     106         609 :       return L? lg(L)-1: 0;
     107             :     }
     108          14 :     case t_REAL: return signe(x)? lg(x)-2: 0;
     109          11 :     case t_STR:  return strlen( GSTR(x) );
     110          91 :     case t_VECSMALL: return lg(x)-1;
     111             :   }
     112     5144776 :   return lg(x) - lontyp[tx];
     113             : }
     114             : 
     115             : long
     116        3878 : gtranslength(GEN x)
     117             : {
     118        3878 :   switch(typ(x))
     119             :   {
     120           0 :     case t_VEC: case t_COL:
     121           0 :       return lg(x)-1;
     122        3878 :     case t_MAT:
     123        3878 :       return lg(x)==1 ? 0: nbrows(x);
     124           0 :     default:
     125           0 :       pari_err_TYPE("trans",x);
     126             :       return 0; /* LCOV_EXCL_LINE */
     127             :   }
     128             : }
     129             : 
     130             : GEN
     131        1862 : matsize(GEN x)
     132             : {
     133        1862 :   long L = lg(x) - 1;
     134        1862 :   switch(typ(x))
     135             :   {
     136           7 :     case t_VEC: return mkvec2s(1, L);
     137           7 :     case t_COL: return mkvec2s(L, 1);
     138        1841 :     case t_MAT: return mkvec2s(L? nbrows(x): 0, L);
     139             :   }
     140           7 :   pari_err_TYPE("matsize",x);
     141             :   return NULL; /* LCOV_EXCL_LINE */
     142             : }
     143             : 
     144             : /*******************************************************************/
     145             : /*                                                                 */
     146             : /*                 CONVERSION GEN --> long                         */
     147             : /*                                                                 */
     148             : /*******************************************************************/
     149             : 
     150             : long
     151          77 : gtolong(GEN x)
     152             : {
     153          77 :   switch(typ(x))
     154             :   {
     155          42 :     case t_INT:
     156          42 :       return itos(x);
     157           7 :     case t_REAL:
     158           7 :       return (long)(rtodbl(x) + 0.5);
     159           7 :     case t_FRAC:
     160           7 :     { pari_sp av = avma; return gc_long(av, itos(ground(x))); }
     161           7 :     case t_COMPLEX:
     162           7 :       if (gequal0(gel(x,2))) return gtolong(gel(x,1)); break;
     163           7 :     case t_QUAD:
     164           7 :       if (gequal0(gel(x,3))) return gtolong(gel(x,2)); break;
     165             :   }
     166           7 :   pari_err_TYPE("gtolong",x);
     167             :   return 0; /* LCOV_EXCL_LINE */
     168             : }
     169             : 
     170             : /*******************************************************************/
     171             : /*                                                                 */
     172             : /*                         COMPARISONS                             */
     173             : /*                                                                 */
     174             : /*******************************************************************/
     175             : static void
     176         189 : chk_true_err()
     177             : {
     178         189 :   GEN E = pari_err_last();
     179         189 :   switch(err_get_num(E))
     180             :   {
     181           0 :     case e_STACK: case e_MEM: case e_ALARM:
     182           0 :       pari_err(0, E); /* rethrow */
     183             :   }
     184         189 : }
     185             : /* x - y == 0 or undefined */
     186             : static int
     187     3198025 : gequal_try(GEN x, GEN y)
     188             : {
     189             :   int i;
     190     3198025 :   pari_CATCH(CATCH_ALL) { chk_true_err(); return 0; }
     191     3198025 :   pari_TRY { i = gequal0(gadd(x, gneg_i(y))); } pari_ENDCATCH;
     192     3197843 :   return i;
     193             : }
     194             : /* x + y == 0 or undefined */
     195             : static int
     196          28 : gmequal_try(GEN x, GEN y)
     197             : {
     198             :   int i;
     199          28 :   pari_CATCH(CATCH_ALL) { chk_true_err(); return 0; }
     200          28 :   pari_TRY { i = gequal0(gadd(x, y)); } pari_ENDCATCH;
     201          21 :   return i;
     202             : }
     203             : 
     204             : int
     205   455364330 : isexactzero(GEN g)
     206             : {
     207             :   long i, lx;
     208   455364330 :   switch (typ(g))
     209             :   {
     210   397049486 :     case t_INT:
     211   397049486 :       return !signe(g);
     212     1242135 :     case t_INTMOD:
     213     1242135 :       return !signe(gel(g,2));
     214    14863240 :     case t_COMPLEX:
     215    14863240 :       return isexactzero(gel(g,1)) && isexactzero(gel(g,2));
     216     8107278 :     case t_FFELT:
     217     8107278 :       return FF_equal0(g);
     218         511 :     case t_QUAD:
     219         511 :       return isexactzero(gel(g,2)) && isexactzero(gel(g,3));
     220      283813 :     case t_POLMOD:
     221      283813 :       return isexactzero(gel(g,2));
     222    12373008 :     case t_POL:
     223    12373008 :       lx = lg(g); /* cater for Mod(0,2)*x^0 */
     224    12373008 :       return lx == 2 || (lx == 3 && isexactzero(gel(g,2)));
     225      462459 :     case t_RFRAC:
     226      462459 :       return isexactzero(gel(g,1)); /* may occur: Mod(0,2)/x */
     227       43435 :     case t_VEC: case t_COL: case t_MAT:
     228       43820 :       for (i=lg(g)-1; i; i--)
     229       43673 :         if (!isexactzero(gel(g,i))) return 0;
     230         147 :       return 1;
     231             :   }
     232    20938965 :   return 0;
     233             : }
     234             : GEN
     235    63715285 : gisexactzero(GEN g)
     236             : {
     237             :   long i, lx;
     238             :   GEN a, b;
     239    63715285 :   switch (typ(g))
     240             :   {
     241    26743481 :     case t_INT:
     242    26743481 :       return !signe(g)? g: NULL;
     243     7065821 :     case t_INTMOD:
     244     7065821 :       return !signe(gel(g,2))? g: NULL;
     245        2632 :     case t_COMPLEX:
     246        2632 :       a = gisexactzero(gel(g,1)); if (!a) return NULL;
     247         616 :       b = gisexactzero(gel(g,2)); if (!b) return NULL;
     248           0 :       return ggcd(a,b);
     249       20608 :     case t_FFELT:
     250       20608 :       return FF_equal0(g)? g: NULL;
     251         518 :     case t_QUAD:
     252         518 :       a = gisexactzero(gel(g,2)); if (!a) return NULL;
     253          77 :       b = gisexactzero(gel(g,3)); if (!b) return NULL;
     254           7 :       return ggcd(a,b);
     255       17032 :     case t_POLMOD:
     256       17032 :       return gisexactzero(gel(g,2));
     257    28072426 :     case t_POL:
     258    28072426 :       lx = lg(g); /* cater for Mod(0,2)*x^0 */
     259    28072426 :       if (lx == 2) return gen_0;
     260    22750737 :       if (lx == 3) return gisexactzero(gel(g,2));
     261    19028727 :       return NULL;
     262     1190293 :     case t_RFRAC:
     263     1190293 :       return gisexactzero(gel(g,1)); /* may occur: Mod(0,2)/x */
     264           0 :     case t_VEC: case t_COL: case t_MAT:
     265           0 :       a = gen_0;
     266           0 :       for (i=lg(g)-1; i; i--)
     267             :       {
     268           0 :         b = gisexactzero(gel(g,i));
     269           0 :         if (!b) return NULL;
     270           0 :         a = ggcd(a, b);
     271             :       }
     272           0 :       return a;
     273             :   }
     274      602474 :   return NULL;
     275             : }
     276             : 
     277             : int
     278   555321210 : isrationalzero(GEN g)
     279             : {
     280             :   long i;
     281   555321210 :   switch (typ(g))
     282             :   {
     283   347074286 :     case t_INT:
     284   347074286 :       return !signe(g);
     285    39373284 :     case t_COMPLEX:
     286    39373284 :       return isintzero(gel(g,1)) && isintzero(gel(g,2));
     287        1428 :     case t_QUAD:
     288        1428 :       return isintzero(gel(g,2)) && isintzero(gel(g,3));
     289      446873 :     case t_POLMOD:
     290      446873 :       return isrationalzero(gel(g,2));
     291    23685344 :     case t_POL: return lg(g) == 2;
     292         133 :     case t_VEC: case t_COL: case t_MAT:
     293         448 :       for (i=lg(g)-1; i; i--)
     294         315 :         if (!isrationalzero(gel(g,i))) return 0;
     295         133 :       return 1;
     296             :   }
     297   144739862 :   return 0;
     298             : }
     299             : 
     300             : int
     301  2292626456 : gequal0(GEN x)
     302             : {
     303  2292626456 :   switch(typ(x))
     304             :   {
     305  2128489843 :     case t_INT: case t_REAL: case t_POL: case t_SER:
     306  2128489843 :       return !signe(x);
     307             : 
     308     7846378 :     case t_INTMOD:
     309     7846378 :       return !signe(gel(x,2));
     310             : 
     311      650970 :     case t_FFELT:
     312      650970 :       return FF_equal0(x);
     313             : 
     314   105341550 :     case t_COMPLEX:
     315             :      /* is 0 iff norm(x) would be 0 (can happen with Re(x) and Im(x) != 0
     316             :       * only if Re(x) and Im(x) are of type t_REAL). See mp.c:addrr().
     317             :       */
     318   105341550 :       if (gequal0(gel(x,1)))
     319             :       {
     320     8112561 :         if (gequal0(gel(x,2))) return 1;
     321     7728968 :         if (typ(gel(x,1))!=t_REAL || typ(gel(x,2))!=t_REAL) return 0;
     322      271211 :         return (expo(gel(x,1))>=expo(gel(x,2)));
     323             :       }
     324    97229180 :       if (gequal0(gel(x,2)))
     325             :       {
     326     1660982 :         if (typ(gel(x,1))!=t_REAL || typ(gel(x,2))!=t_REAL) return 0;
     327     1576594 :         return (expo(gel(x,2))>=expo(gel(x,1)));
     328             :       }
     329    95572278 :       return 0;
     330             : 
     331     2046636 :     case t_PADIC:
     332     2046636 :       return !signe(gel(x,4));
     333             : 
     334        1806 :     case t_QUAD:
     335        1806 :       return gequal0(gel(x,2)) && gequal0(gel(x,3));
     336             : 
     337     8507875 :     case t_POLMOD:
     338     8507875 :       return gequal0(gel(x,2));
     339             : 
     340     6104264 :     case t_RFRAC:
     341     6104264 :       return gequal0(gel(x,1));
     342             : 
     343     9919814 :     case t_VEC: case t_COL: case t_MAT:
     344             :     {
     345             :       long i;
     346    22873231 :       for (i=lg(x)-1; i; i--)
     347    19091154 :         if (!gequal0(gel(x,i))) return 0;
     348     3782077 :       return 1;
     349             :     }
     350             :   }
     351    23717320 :   return 0;
     352             : }
     353             : 
     354             : /* x a t_POL or t_SER, return 1 if test(coeff(X,d)) is true and
     355             :  * coeff(X,i) = 0 for all i != d. Return 0 (false) otherwise */
     356             : static int
     357    15238306 : is_monomial_test(GEN x, long d, int(*test)(GEN))
     358             : {
     359    15238306 :   long i, l = lg(x);
     360    15238306 :   if (typ(x) == t_SER)
     361             :   { /* "0" * x^v * (1+O(x)) ?  v <= 0 or null ring */
     362         602 :     if (l == 3 && isexactzero(gel(x,2))) return d >= 2 || test(gel(x,2));
     363         553 :     if (d < 2) return 0; /* v > 0 */
     364             :   }
     365    15238075 :   if (d >= l)
     366             :   {
     367       60086 :     if (typ(x) == t_POL) return 0; /* l = 2 */
     368             :     /* t_SER, v = 2-d <= 0 */
     369          56 :     if (!signe(x)) return 1;
     370             :   }
     371    15177989 :   else if (!test(gel(x,d))) return 0;
     372     7357684 :   for (i = 2; i < l; i++) /* 2 <= d < l */
     373     4977336 :     if (i != d && !gequal0(gel(x,i))) return 0;
     374     2380348 :   return 1;
     375             : }
     376             : static int
     377      294364 : col_test(GEN x, int(*test)(GEN))
     378             : {
     379      294364 :   long i, l = lg(x);
     380      294364 :   if (l == 1 || !test(gel(x,1))) return 0;
     381       13167 :   for (i = 2; i < l; i++)
     382       11186 :     if (!gequal0(gel(x,i))) return 0;
     383        1981 :   return 1;
     384             : }
     385             : static int
     386       16310 : mat_test(GEN x, int(*test)(GEN))
     387             : {
     388       16310 :   long i, j, l = lg(x);
     389       16310 :   if (l == 1) return 1;
     390       16296 :   if (l != lgcols(x)) return 0;
     391       51779 :   for (i = 1; i < l; i++)
     392      136941 :     for (j = 1; j < l; j++)
     393      101458 :       if (i == j) {
     394       35567 :         if (!test(gcoeff(x,i,i))) return 0;
     395             :       } else {
     396       65891 :         if (!gequal0(gcoeff(x,i,j))) return 0;
     397             :       }
     398       16212 :   return 1;
     399             : }
     400             : 
     401             : /* returns 1 whenever x = 1, and 0 otherwise */
     402             : int
     403   328384279 : gequal1(GEN x)
     404             : {
     405   328384279 :   switch(typ(x))
     406             :   {
     407   308433573 :     case t_INT:
     408   308433573 :       return equali1(x);
     409             : 
     410       71050 :     case t_REAL:
     411             :     {
     412       71050 :       long s = signe(x);
     413       71050 :       if (!s) return expo(x) >= 0;
     414       70952 :       return s > 0 ? absrnz_equal1(x): 0;
     415             :     }
     416      632454 :     case t_INTMOD:
     417      632454 :       return is_pm1(gel(x,2)) || is_pm1(gel(x,1));
     418      393547 :     case t_POLMOD:
     419      393547 :       return !degpol(gel(x,1)) || gequal1(gel(x,2));
     420             : 
     421       16891 :     case t_FFELT:
     422       16891 :       return FF_equal1(x);
     423             : 
     424     1699608 :     case t_FRAC:
     425     1699608 :       return 0;
     426             : 
     427       23994 :     case t_COMPLEX:
     428       23994 :       return gequal1(gel(x,1)) && gequal0(gel(x,2));
     429             : 
     430      166225 :     case t_PADIC:
     431      166225 :       if (!signe(gel(x,4))) return valp(x) <= 0;
     432      166183 :       return valp(x) == 0 && gequal1(gel(x,4));
     433             : 
     434          42 :     case t_QUAD:
     435          42 :       return gequal1(gel(x,2)) && gequal0(gel(x,3));
     436             : 
     437    15237634 :     case t_POL: return is_monomial_test(x, 2, &gequal1);
     438         476 :     case t_SER: return is_monomial_test(x, 2 - valser(x), &gequal1);
     439             : 
     440     1028431 :     case t_RFRAC: return gequal(gel(x,1), gel(x,2));
     441      294315 :     case t_COL: return col_test(x, &gequal1);
     442       16198 :     case t_MAT: return mat_test(x, &gequal1);
     443             :   }
     444      369841 :   return 0;
     445             : }
     446             : 
     447             : /* returns 1 whenever the x = -1, 0 otherwise */
     448             : int
     449    74232958 : gequalm1(GEN x)
     450             : {
     451             :   pari_sp av;
     452             :   GEN t;
     453             : 
     454    74232958 :   switch(typ(x))
     455             :   {
     456    74224083 :     case t_INT:
     457    74224083 :       return equalim1(x);
     458             : 
     459        1484 :     case t_REAL:
     460             :     {
     461        1484 :       long s = signe(x);
     462        1484 :       if (!s) return expo(x) >= 0;
     463        1477 :       return s < 0 ? absrnz_equal1(x): 0;
     464             :     }
     465        4640 :     case t_INTMOD:
     466        4640 :       av = avma; return gc_bool(av, equalii(addui(1,gel(x,2)), gel(x,1)));
     467             : 
     468         154 :     case t_FRAC:
     469         154 :       return 0;
     470             : 
     471          42 :     case t_FFELT:
     472          42 :       return FF_equalm1(x);
     473             : 
     474        2044 :     case t_COMPLEX:
     475        2044 :       return gequalm1(gel(x,1)) && gequal0(gel(x,2));
     476             : 
     477           7 :     case t_QUAD:
     478           7 :       return gequalm1(gel(x,2)) && gequal0(gel(x,3));
     479             : 
     480          49 :     case t_PADIC:
     481          49 :       t = gel(x,4); if (!signe(t)) return valp(x) <= 0;
     482          21 :       av = avma; return gc_bool(av, !valp(x) && equalii(addui(1,t), gel(x,3)));
     483             : 
     484          56 :     case t_POLMOD:
     485          56 :       return !degpol(gel(x,1)) || gequalm1(gel(x,2));
     486             : 
     487          70 :     case t_POL: return is_monomial_test(x, 2, &gequalm1);
     488         126 :     case t_SER: return is_monomial_test(x, 2 - valser(x), &gequalm1);
     489             : 
     490          28 :     case t_RFRAC:
     491          28 :       av = avma; return gc_bool(av, gmequal_try(gel(x,1), gel(x,2)));
     492          49 :     case t_COL: return col_test(x, &gequalm1);
     493         112 :     case t_MAT: return mat_test(x, &gequalm1);
     494             :   }
     495          14 :   return 0;
     496             : }
     497             : 
     498             : int
     499     1469454 : gequalX(GEN x) { return typ(x) == t_POL && lg(x) == 4
     500     2953768 :                       && isintzero(gel(x,2)) && isint1(gel(x,3)); }
     501             : 
     502             : static int
     503         672 : cmp_str(const char *x, const char *y)
     504             : {
     505         672 :   int f = strcmp(x, y);
     506             :   return f > 0? 1
     507         672 :               : f? -1: 0;
     508             : }
     509             : 
     510             : static int
     511    39216142 : cmp_universal_rec(GEN x, GEN y, long i0)
     512             : {
     513    39216142 :   long i, lx = lg(x), ly = lg(y);
     514    39216142 :   if (lx < ly) return -1;
     515    39213367 :   if (lx > ly) return 1;
     516    69274401 :   for (i = i0; i < lx; i++)
     517             :   {
     518    60359648 :     int f = cmp_universal(gel(x,i), gel(y,i));
     519    60359648 :     if (f) return f;
     520             :   }
     521     8914753 :   return 0;
     522             : }
     523             : /* Universal "meaningless" comparison function. Transitive, returns 0 iff
     524             :  * gidentical(x,y) */
     525             : int
     526    84946542 : cmp_universal(GEN x, GEN y)
     527             : {
     528    84946542 :   long lx, ly, i, tx = typ(x), ty = typ(y);
     529             : 
     530    84946542 :   if (tx < ty) return -1;
     531    84573734 :   if (ty < tx) return 1;
     532    84045036 :   switch(tx)
     533             :   {
     534    43893161 :     case t_INT: return cmpii(x,y);
     535         651 :     case t_STR: return cmp_str(GSTR(x),GSTR(y));
     536      934997 :     case t_REAL:
     537             :     case t_VECSMALL:
     538      934997 :       lx = lg(x);
     539      934997 :       ly = lg(y);
     540      934997 :       if (lx < ly) return -1;
     541      886725 :       if (lx > ly) return 1;
     542     3586440 :       for (i = 1; i < lx; i++)
     543             :       {
     544     3478808 :         if (x[i] < y[i]) return -1;
     545     3109879 :         if (x[i] > y[i]) return 1;
     546             :       }
     547      107632 :       return 0;
     548             : 
     549      771675 :     case t_POL:
     550             :     {
     551      771675 :       long X = x[1] & (VARNBITS|SIGNBITS);
     552      771675 :       long Y = y[1] & (VARNBITS|SIGNBITS);
     553      771675 :       if (X < Y) return -1;
     554      771654 :       if (X > Y) return 1;
     555      771598 :       return cmp_universal_rec(x, y, 2);
     556             :     }
     557      881076 :     case t_SER:
     558             :     case t_FFELT:
     559             :     case t_CLOSURE:
     560      881076 :       if (x[1] < y[1]) return -1;
     561      881069 :       if (x[1] > y[1]) return 1;
     562      881062 :       return cmp_universal_rec(x, y, 2);
     563             : 
     564          35 :     case t_LIST:
     565             :       {
     566          35 :         long tx = list_typ(x), ty = list_typ(y);
     567             :         GEN vx, vy;
     568             :         pari_sp av;
     569          35 :         if (tx < ty) return -1;
     570          35 :         if (tx > ty) return 1;
     571          35 :         vx = list_data(x);
     572          35 :         vy = list_data(y);
     573          35 :         if (!vx) return vy? -1: 0;
     574          35 :         if (!vy) return 1;
     575          35 :         av = avma;
     576          35 :         if (tx == t_LIST_MAP)
     577             :         {
     578          14 :           vx = maptomat_shallow(x);
     579          14 :           vy = maptomat_shallow(y);
     580             :         }
     581          35 :         return gc_int(av, cmp_universal_rec(vx, vy, 1));
     582             :       }
     583    37563441 :     default:
     584    37563441 :       return cmp_universal_rec(x, y, lontyp[tx]);
     585             :   }
     586             : }
     587             : 
     588             : static int
     589     4588920 : cmpfrac(GEN x, GEN y)
     590             : {
     591     4588920 :   pari_sp av = avma;
     592     4588920 :   GEN a = gel(x,1), b = gel(x,2);
     593     4588920 :   GEN c = gel(y,1), d = gel(y,2);
     594     4588920 :   return gc_bool(av, cmpii(mulii(a, d), mulii(b, c)));
     595             : }
     596             : static int
     597      455541 : cmpifrac(GEN a, GEN y)
     598             : {
     599      455541 :   pari_sp av = avma;
     600      455541 :   GEN c = gel(y,1), d = gel(y,2);
     601      455541 :   return gc_int(av, cmpii(mulii(a, d), c));
     602             : }
     603             : static int
     604       81189 : cmprfrac(GEN a, GEN y)
     605             : {
     606       81189 :   pari_sp av = avma;
     607       81189 :   GEN c = gel(y,1), d = gel(y,2);
     608       81189 :   return gc_int(av, cmpri(mulri(a, d), c));
     609             : }
     610             : static int
     611         161 : cmpgen(GEN x, GEN y)
     612             : {
     613         161 :   pari_sp av = avma;
     614         161 :   return gc_int(av, gsigne(gsub(x,y)));
     615             : }
     616             : 
     617             : /* returns the sign of x - y when it makes sense. 0 otherwise */
     618             : int
     619   290887277 : gcmp(GEN x, GEN y)
     620             : {
     621   290887277 :   long tx = typ(x), ty = typ(y);
     622             : 
     623   290887277 :   if (tx == ty) /* generic case */
     624   281376618 :     switch(tx)
     625             :     {
     626   153353542 :       case t_INT:  return cmpii(x, y);
     627   123339948 :       case t_REAL: return cmprr(x, y);
     628     4588920 :       case t_FRAC: return cmpfrac(x, y);
     629          70 :       case t_QUAD: return cmpgen(x, y);
     630          21 :       case t_STR:  return cmp_str(GSTR(x), GSTR(y));
     631      104551 :       case t_INFINITY:
     632             :       {
     633      104551 :         long sx = inf_get_sign(x), sy = inf_get_sign(y);
     634      104551 :         if (sx < sy) return -1;
     635          42 :         if (sx > sy) return 1;
     636          14 :         return 0;
     637             :       }
     638             :     }
     639     9500225 :   if (ty == t_INFINITY) return -inf_get_sign(y);
     640     8555097 :   switch(tx)
     641             :   {
     642     8134044 :     case t_INT:
     643             :       switch(ty)
     644             :       {
     645     7821268 :         case t_REAL: return cmpir(x, y);
     646      312761 :         case t_FRAC: return cmpifrac(x, y);
     647           7 :         case t_QUAD: return cmpgen(x, y);
     648             :       }
     649           8 :       break;
     650      215160 :     case t_REAL:
     651             :       switch(ty)
     652             :       {
     653      178292 :         case t_INT:  return cmpri(x, y);
     654       36847 :         case t_FRAC: return cmprfrac(x, y);
     655          14 :         case t_QUAD: return cmpgen(x, y);
     656             :       }
     657           7 :       break;
     658      187136 :     case t_FRAC:
     659             :       switch(ty)
     660             :       {
     661      142780 :         case t_INT:  return -cmpifrac(y, x);
     662       44342 :         case t_REAL: return -cmprfrac(y, x);
     663           7 :         case t_QUAD: return cmpgen(x, y);
     664             :       }
     665           7 :       break;
     666          63 :     case t_QUAD:
     667          63 :       return cmpgen(x, y);
     668       31646 :     case t_INFINITY: return inf_get_sign(x);
     669             :   }
     670          24 :   pari_err_TYPE2("comparison",x,y);
     671             :   return 0;/*LCOV_EXCL_LINE*/
     672             : }
     673             : 
     674             : int
     675      780700 : gcmpsg(long s, GEN y)
     676             : {
     677      780700 :   switch(typ(y))
     678             :   {
     679       12236 :     case t_INT:  return cmpsi(s,y);
     680      763291 :     case t_REAL: return cmpsr(s,y);
     681        5173 :     case t_FRAC: {
     682        5173 :       pari_sp av = avma;
     683        5173 :       return gc_int(av, cmpii(mulsi(s,gel(y,2)), gel(y,1)));
     684             :     }
     685           0 :     case t_QUAD: {
     686           0 :       pari_sp av = avma;
     687           0 :       return gc_int(av, gsigne(gsubsg(s, y)));
     688             :     }
     689           0 :     case t_INFINITY: return -inf_get_sign(y);
     690             :   }
     691           0 :   pari_err_TYPE2("comparison",stoi(s),y);
     692             :   return 0; /* LCOV_EXCL_LINE */
     693             : }
     694             : 
     695             : static long
     696     3232583 : roughtype(GEN x)
     697             : {
     698     3232583 :   switch(typ(x))
     699             :   {
     700        2114 :     case t_MAT: return t_MAT;
     701      742362 :     case t_VEC: case t_COL: return t_VEC;
     702     1613554 :     case t_VECSMALL: return t_VECSMALL;
     703      874553 :     default: return t_INT;
     704             :   }
     705             : }
     706             : 
     707             : static int lexcmpsg(long x, GEN y);
     708          42 : static int lexcmpgs(GEN x, long y) { return -lexcmpsg(y,x); }
     709             : /* lexcmp(stoi(x),y), y t_VEC/t_COL/t_MAT */
     710             : static int
     711          21 : lexcmp_s_matvec(long x, GEN y)
     712             : {
     713             :   int fl;
     714          21 :   if (lg(y)==1) return 1;
     715          14 :   fl = lexcmpsg(x,gel(y,1));
     716          14 :   if (fl) return fl;
     717           7 :   return -1;
     718             : }
     719             : /* x a scalar, y a t_VEC/t_COL/t_MAT */
     720             : static int
     721         357 : lexcmp_scal_matvec(GEN x, GEN y)
     722             : {
     723             :   int fl;
     724         357 :   if (lg(y)==1) return 1;
     725         357 :   fl = lexcmp(x,gel(y,1));
     726         357 :   if (fl) return fl;
     727           7 :   return -1;
     728             : }
     729             : /* x a scalar, y a t_VECSMALL */
     730             : static int
     731          42 : lexcmp_scal_vecsmall(GEN x, GEN y)
     732             : {
     733             :   int fl;
     734          42 :   if (lg(y)==1) return 1;
     735          42 :   fl = lexcmpgs(x, y[1]);
     736          42 :   if (fl) return fl;
     737           0 :   return -1;
     738             : }
     739             : 
     740             : /* tx = ty = t_MAT, or x and y are both vect_t */
     741             : static int
     742      371973 : lexcmp_similar(GEN x, GEN y)
     743             : {
     744      371973 :   long i, lx = lg(x), ly = lg(y), l = minss(lx,ly);
     745      456824 :   for (i=1; i<l; i++)
     746             :   {
     747      425199 :     int fl = lexcmp(gel(x,i),gel(y,i));
     748      425198 :     if (fl) return fl;
     749             :   }
     750       31625 :   if (lx == ly) return 0;
     751          34 :   return (lx < ly)? -1 : 1;
     752             : }
     753             : /* x a t_VECSMALL, y a t_VEC/t_COL ~ lexcmp_similar */
     754             : static int
     755         154 : lexcmp_vecsmall_vec(GEN x, GEN y)
     756             : {
     757         154 :   long i, lx = lg(x), ly = lg(y), l = minss(lx,ly);
     758         343 :   for (i=1; i<l; i++)
     759             :   {
     760         287 :     int fl = lexcmpsg(x[i], gel(y,i));
     761         287 :     if (fl) return fl;
     762             :   }
     763          56 :   if (lx == ly) return 0;
     764          21 :   return (lx < ly)? -1 : 1;
     765             : }
     766             : 
     767             : /* x t_VEC/t_COL, y t_MAT */
     768             : static int
     769          98 : lexcmp_vec_mat(GEN x, GEN y)
     770             : {
     771             :   int fl;
     772          98 :   if (lg(x)==1) return -1;
     773          98 :   if (lg(y)==1) return 1;
     774          98 :   fl = lexcmp_similar(x,gel(y,1));
     775          98 :   if (fl) return fl;
     776           7 :   return -1;
     777             : }
     778             : /* x t_VECSMALl, y t_MAT ~ lexcmp_vec_mat */
     779             : static int
     780          42 : lexcmp_vecsmall_mat(GEN x, GEN y)
     781             : {
     782             :   int fl;
     783          42 :   if (lg(x)==1) return -1;
     784          42 :   if (lg(y)==1) return 1;
     785          42 :   fl = lexcmp_vecsmall_vec(x, gel(y,1));
     786          42 :   if (fl) return fl;
     787           0 :   return -1;
     788             : }
     789             : 
     790             : /* x a t_VECSMALL, not y */
     791             : static int
     792         196 : lexcmp_vecsmall_other(GEN x, GEN y, long ty)
     793             : {
     794         196 :   switch(ty)
     795             :   {
     796          42 :     case t_MAT: return lexcmp_vecsmall_mat(x, y);
     797         112 :     case t_VEC: return lexcmp_vecsmall_vec(x, y);
     798          42 :     default: return -lexcmp_scal_vecsmall(y, x); /*y scalar*/
     799             :   }
     800             : }
     801             : 
     802             : /* lexcmp(stoi(s), y) */
     803             : static int
     804         343 : lexcmpsg(long x, GEN y)
     805             : {
     806         343 :   switch(roughtype(y))
     807             :   {
     808          21 :     case t_MAT:
     809             :     case t_VEC:
     810          21 :       return lexcmp_s_matvec(x,y);
     811          14 :     case t_VECSMALL: /* ~ lexcmp_scal_matvec */
     812          14 :       if (lg(y)==1) return 1;
     813           7 :       return (x > y[1])? 1: -1;
     814         308 :     default: return gcmpsg(x,y);
     815             :   }
     816             : }
     817             : 
     818             : /* as gcmp for vector/matrices, using lexicographic ordering on components */
     819             : static int
     820     1616120 : lexcmp_i(GEN x, GEN y)
     821             : {
     822     1616120 :   const long tx = roughtype(x), ty = roughtype(y);
     823     1616121 :   if (tx == ty)
     824     1615470 :     switch(tx)
     825             :     {
     826      371874 :       case t_MAT:
     827      371874 :       case t_VEC: return lexcmp_similar(x,y);
     828      806672 :       case t_VECSMALL: return vecsmall_lexcmp(x,y);
     829      436924 :       default: return gcmp(x,y);
     830             :     }
     831         651 :   if (tx == t_VECSMALL) return  lexcmp_vecsmall_other(x,y,ty);
     832         518 :   if (ty == t_VECSMALL) return -lexcmp_vecsmall_other(y,x,tx);
     833             : 
     834         455 :   if (tx == t_INT) return  lexcmp_scal_matvec(x,y); /*scalar*/
     835         203 :   if (ty == t_INT) return -lexcmp_scal_matvec(y,x);
     836             : 
     837          98 :   if (ty==t_MAT) return  lexcmp_vec_mat(x,y);
     838          42 :   return -lexcmp_vec_mat(y,x); /*tx==t_MAT*/
     839             : }
     840             : int
     841     1616120 : lexcmp(GEN x, GEN y)
     842             : {
     843     1616120 :   pari_sp av = avma;
     844     1616120 :   if (typ(x) == t_COMPLEX)
     845             :   {
     846         875 :     x = mkvec2(gel(x,1), gel(x,2));
     847         875 :     if (typ(y) == t_COMPLEX) y = mkvec2(gel(y,1), gel(y,2));
     848          49 :     else y = mkvec2(y, gen_0);
     849             :   }
     850     1615245 :   else if (typ(y) == t_COMPLEX)
     851             :   {
     852          63 :     x = mkvec2(x, gen_0);
     853          63 :     y = mkvec2(gel(y,1), gel(y,2));
     854             :   }
     855     1616120 :   return gc_int(av, lexcmp_i(x, y));
     856             : }
     857             : 
     858             : /*****************************************************************/
     859             : /*                                                               */
     860             : /*                          EQUALITY                             */
     861             : /*                returns 1 if x == y, 0 otherwise               */
     862             : /*                                                               */
     863             : /*****************************************************************/
     864             : /* x,y t_POL */
     865             : static int
     866     3432210 : polidentical(GEN x, GEN y)
     867             : {
     868             :   long lx;
     869     3432210 :   if (x[1] != y[1]) return 0;
     870     3432112 :   lx = lg(x); if (lg(y) != lg(x)) return 0;
     871    14948876 :   for (lx--; lx >= 2; lx--) if (!gidentical(gel(x,lx), gel(y,lx))) return 0;
     872     3432021 :   return 1;
     873             : }
     874             : /* x,y t_SER */
     875             : static int
     876          14 : seridentical(GEN x, GEN y) { return polidentical(x,y); }
     877             : /* typ(x) = typ(y) = t_VEC/COL/MAT */
     878             : static int
     879     5264835 : vecidentical(GEN x, GEN y)
     880             : {
     881             :   long i;
     882     5264835 :   if ((x[0] ^ y[0]) & (TYPBITS|LGBITS)) return 0;
     883    16813713 :   for (i = lg(x)-1; i; i--)
     884    12843297 :     if (! gidentical(gel(x,i),gel(y,i)) ) return 0;
     885     3970416 :   return 1;
     886             : }
     887             : static int
     888        1547 : identicalrr(GEN x, GEN y)
     889             : {
     890        1547 :   long i, lx = lg(x);
     891        1547 :   if (lg(y) != lx) return 0;
     892        1547 :   if (x[1] != y[1]) return 0;
     893        5465 :   i=2; while (i<lx && x[i]==y[i]) i++;
     894        1540 :   return (i == lx);
     895             : }
     896             : 
     897             : static int
     898          70 : closure_identical(GEN x, GEN y)
     899             : {
     900          70 :   if (lg(x)!=lg(y) || x[1]!=y[1]) return 0;
     901          56 :   if (!gidentical(gel(x,2),gel(y,2)) || !gidentical(gel(x,3),gel(y,3))
     902          56 :    || !gidentical(gel(x,4),gel(y,4))) return 0;
     903          42 :   if (lg(x)<8) return 1;
     904           0 :   return gidentical(gel(x,7),gel(y,7));
     905             : }
     906             : 
     907             : static int
     908         343 : list_cmp(GEN x, GEN y, int cmp(GEN x, GEN y))
     909             : {
     910         343 :   int t = list_typ(x);
     911             :   GEN vx, vy;
     912             :   long lvx, lvy;
     913         343 :   if (list_typ(y)!=t) return 0;
     914         343 :   vx = list_data(x);
     915         343 :   vy = list_data(y);
     916         343 :   lvx = vx ? lg(vx): 1;
     917         343 :   lvy = vy ? lg(vy): 1;
     918         343 :   if (lvx==1 && lvy==1) return 1;
     919         329 :   if (lvx != lvy) return 0;
     920         301 :   switch (t)
     921             :   {
     922         280 :   case t_LIST_MAP:
     923             :     {
     924         280 :       pari_sp av = avma;
     925         280 :       GEN mx  = maptomat_shallow(x), my = maptomat_shallow(y);
     926         280 :       int ret = gidentical(gel(mx, 1), gel(my, 1)) && cmp(gel(mx, 2), gel(my, 2));
     927         280 :       return gc_bool(av, ret);
     928             :     }
     929          21 :   default:
     930          21 :     return cmp(vx, vy);
     931             :   }
     932             : }
     933             : 
     934             : int
     935    62038119 : gidentical(GEN x, GEN y)
     936             : {
     937             :   long tx;
     938             : 
     939    62038119 :   if (x == y) return 1;
     940    58333821 :   tx = typ(x); if (typ(y) != tx) return 0;
     941    58099976 :   switch(tx)
     942             :   {
     943    19082922 :     case t_INT:
     944    19082922 :       return equalii(x,y);
     945             : 
     946        1547 :     case t_REAL:
     947        1547 :       return identicalrr(x,y);
     948             : 
     949      794026 :     case t_FRAC: case t_INTMOD:
     950      794026 :       return equalii(gel(x,2), gel(y,2)) && equalii(gel(x,1), gel(y,1));
     951             : 
     952         343 :     case t_COMPLEX:
     953         343 :       return gidentical(gel(x,2),gel(y,2)) && gidentical(gel(x,1),gel(y,1));
     954          14 :     case t_PADIC:
     955          14 :       return valp(x) == valp(y)
     956          14 :         && equalii(gel(x,2),gel(y,2))
     957          14 :         && equalii(gel(x,3),gel(y,3))
     958          28 :         && equalii(gel(x,4),gel(y,4));
     959        3913 :     case t_POLMOD:
     960        3913 :       return gidentical(gel(x,2),gel(y,2)) && polidentical(gel(x,1),gel(y,1));
     961     3432175 :     case t_POL:
     962     3432175 :       return polidentical(x,y);
     963          14 :     case t_SER:
     964          14 :       return seridentical(x,y);
     965        3024 :     case t_FFELT:
     966        3024 :       return FF_equal(x,y);
     967             : 
     968      401624 :     case t_QFB:
     969      401624 :       return equalii(gel(x,1),gel(y,1))
     970      401617 :           && equalii(gel(x,2),gel(y,2))
     971      803241 :           && equalii(gel(x,3),gel(y,3));
     972             : 
     973          14 :     case t_QUAD:
     974          14 :       return ZX_equal(gel(x,1),gel(y,1))
     975           7 :           && gidentical(gel(x,2),gel(y,2))
     976          21 :           && gidentical(gel(x,3),gel(y,3));
     977             : 
     978           7 :     case t_RFRAC:
     979           7 :       return gidentical(gel(x,1),gel(y,1)) && gidentical(gel(x,2),gel(y,2));
     980             : 
     981          70 :     case t_STR:
     982          70 :       return !strcmp(GSTR(x),GSTR(y));
     983     5264835 :     case t_VEC: case t_COL: case t_MAT:
     984     5264835 :       return vecidentical(x,y);
     985    29115238 :     case t_VECSMALL:
     986    29115238 :       return zv_equal(x,y);
     987          28 :     case t_CLOSURE:
     988          28 :       return closure_identical(x,y);
     989         161 :     case t_LIST:
     990         161 :       return list_cmp(x, y, gidentical);
     991          21 :     case t_INFINITY: return gidentical(gel(x,1),gel(y,1));
     992             :   }
     993           0 :   return 0;
     994             : }
     995             : /* x,y t_POL in the same variable */
     996             : static int
     997     7684368 : polequal(GEN x, GEN y)
     998             : {
     999             :   long lx, ly;
    1000             :   /* Can't do that: Mod(0,1)*x^0 == x^0
    1001             :   if (signe(x) != signe(y)) return 0; */
    1002     7684368 :   lx = lg(x); ly = lg(y);
    1003     7684368 :   while (lx > ly) if (!gequal0(gel(x,--lx))) return 0;
    1004     7680819 :   while (ly > lx) if (!gequal0(gel(y,--ly))) return 0;
    1005    30564814 :   for (lx--; lx >= 2; lx--) if (!gequal(gel(x,lx), gel(y,lx))) return 0;
    1006     7608792 :   return 1;
    1007             : }
    1008             : 
    1009             : /* x,y t_SER in the same variable */
    1010             : static int
    1011         420 : serequal(GEN x, GEN y)
    1012             : {
    1013             :   long LX, LY, lx, ly, vx, vy;
    1014         420 :   if (!signe(x) && !signe(y)) return 1;
    1015          56 :   lx = lg(x); vx = valser(x); LX = lx + vx;
    1016          56 :   ly = lg(y); vy = valser(y); LY = ly + vy;
    1017          56 :   if (LX > LY) lx = LY - vx; else ly = LX - vy;
    1018      282877 :   while (lx >= 3 && ly >= 3)
    1019      282821 :     if (!gequal(gel(x,--lx), gel(y,--ly))) return 0;
    1020          56 :   while(--ly >= 2) if (!gequal0(gel(y,ly))) return 0;
    1021          84 :   while(--lx >= 2) if (!gequal0(gel(x,lx))) return 0;
    1022          49 :   return 1;
    1023             : }
    1024             : 
    1025             : /* typ(x) = typ(y) = t_VEC/COL/MAT */
    1026             : static int
    1027     5652569 : vecequal(GEN x, GEN y)
    1028             : {
    1029             :   long i;
    1030     5652569 :   if ((x[0] ^ y[0]) & (TYPBITS|LGBITS)) return 0;
    1031    18780028 :   for (i = lg(x)-1; i; i--)
    1032    16374674 :     if (! gequal(gel(x,i),gel(y,i)) ) return 0;
    1033     2405354 :   return 1;
    1034             : }
    1035             : 
    1036             : int
    1037   234110836 : gequal(GEN x, GEN y)
    1038             : {
    1039             :   pari_sp av;
    1040             :   long tx, ty;
    1041             :   long i;
    1042             : 
    1043   234110836 :   if (x == y) return 1;
    1044   204099992 :   tx = typ(x);
    1045   204099992 :   ty = typ(y);
    1046   204099992 :   if (tx == ty)
    1047   196302957 :     switch(tx)
    1048             :     {
    1049   170187819 :       case t_INT:
    1050   170187819 :         return equalii(x,y);
    1051             : 
    1052       20425 :       case t_REAL:
    1053       20425 :         return equalrr(x,y);
    1054             : 
    1055     6844975 :       case t_FRAC: case t_INTMOD:
    1056     6844975 :         return equalii(gel(x,2), gel(y,2)) && equalii(gel(x,1), gel(y,1));
    1057             : 
    1058        1309 :       case t_COMPLEX:
    1059        1309 :         return gequal(gel(x,2),gel(y,2)) && gequal(gel(x,1),gel(y,1));
    1060         770 :       case t_PADIC:
    1061         770 :         if (!equalii(gel(x,2),gel(y,2))) return 0;
    1062         770 :         av = avma; i = gequal0(gsub(x,y)); set_avma(av);
    1063         770 :         return i;
    1064     3188275 :       case t_POLMOD:
    1065     3188275 :         if (varn(gel(x,1)) != varn(gel(y,1))) break;
    1066     3188254 :         return gequal(gel(x,2),gel(y,2)) && RgX_equal_var(gel(x,1),gel(y,1));
    1067     7692887 :       case t_POL:
    1068     7692887 :         if (varn(x) != varn(y)) break;
    1069     7684368 :         return polequal(x,y);
    1070         420 :       case t_SER:
    1071         420 :         if (varn(x) != varn(y)) break;
    1072         420 :         return serequal(x,y);
    1073             : 
    1074       61250 :       case t_FFELT:
    1075       61250 :         return FF_equal(x,y);
    1076             : 
    1077     1097034 :       case t_QFB:
    1078     1097034 :         return equalii(gel(x,1),gel(y,1))
    1079      247773 :             && equalii(gel(x,2),gel(y,2))
    1080     1344807 :             && equalii(gel(x,3),gel(y,3));
    1081             : 
    1082           7 :       case t_QUAD:
    1083           7 :         return ZX_equal(gel(x,1),gel(y,1))
    1084           0 :             && gequal(gel(x,2),gel(y,2))
    1085           7 :             && gequal(gel(x,3),gel(y,3));
    1086             : 
    1087       73759 :       case t_RFRAC:
    1088             :       {
    1089       73759 :         GEN a = gel(x,1), b = gel(x,2), c = gel(y,1), d = gel(y,2);
    1090       73759 :         if (gequal(b,d)) return gequal(a,c); /* simple case */
    1091           0 :         av = avma;
    1092           0 :         a = simplify_shallow(gmul(a,d));
    1093           0 :         b = simplify_shallow(gmul(b,c));
    1094           0 :         return gc_bool(av, gequal(a,b));
    1095             :       }
    1096             : 
    1097       64939 :       case t_STR:
    1098       64939 :         return !strcmp(GSTR(x),GSTR(y));
    1099     5652569 :       case t_VEC: case t_COL: case t_MAT:
    1100     5652569 :         return vecequal(x,y);
    1101     1416270 :       case t_VECSMALL:
    1102     1416270 :         return zv_equal(x,y);
    1103         182 :       case t_LIST:
    1104         182 :         return list_cmp(x, y, gequal);
    1105          42 :       case t_CLOSURE:
    1106          42 :         return closure_identical(x,y);
    1107          28 :       case t_INFINITY:
    1108          28 :         return gequal(gel(x,1),gel(y,1));
    1109             :     }
    1110     7805572 :   if (is_noncalc_t(tx) || is_noncalc_t(ty)) return 0;
    1111     7805551 :   if (tx == t_INT && !signe(x)) return gequal0(y);
    1112     7802303 :   if (ty == t_INT && !signe(y)) return gequal0(x);
    1113     3198025 :   (void)&av; av = avma; /* emulate volatile */
    1114     3198025 :   return gc_bool(av, gequal_try(x, y));
    1115             : }
    1116             : 
    1117             : int
    1118       43988 : gequalsg(long s, GEN x)
    1119       43988 : { pari_sp av = avma; return gc_bool(av, gequal(stoi(s), x)); }
    1120             : 
    1121             : /* a and b are t_INT, t_FRAC, t_REAL or t_COMPLEX of those. Check whether
    1122             :  * a-b is invertible */
    1123             : int
    1124       53397 : cx_approx_equal(GEN a, GEN b)
    1125             : {
    1126       53397 :   pari_sp av = avma;
    1127             :   GEN d;
    1128       53397 :   if (a == b) return 1;
    1129       24486 :   d = gsub(a,b);
    1130       24486 :   return gc_bool(av, gequal0(d) || (typ(d)==t_COMPLEX && gequal0(cxnorm(d))));
    1131             : }
    1132             : static int
    1133     1748650 : r_approx0(GEN x, long e) { return e - expo(x) > bit_prec(x); }
    1134             : /* x ~ 0 compared to reference y */
    1135             : int
    1136     2477529 : cx_approx0(GEN x, GEN y)
    1137             : {
    1138             :   GEN a, b;
    1139             :   long e;
    1140     2477529 :   switch(typ(x))
    1141             :   {
    1142         469 :     case t_COMPLEX:
    1143         469 :       a = gel(x,1); b = gel(x,2);
    1144         469 :       if (typ(a) != t_REAL)
    1145             :       {
    1146          14 :         if (!gequal0(a)) return 0;
    1147           0 :         a = NULL;
    1148             :       }
    1149         455 :       else if (!signe(a)) a = NULL;
    1150         455 :       if (typ(b) != t_REAL)
    1151             :       {
    1152           0 :         if (!gequal0(b)) return 0;
    1153           0 :         if (!a) return 1;
    1154           0 :         b = NULL;
    1155             :       }
    1156         455 :       else if (!signe(b))
    1157             :       {
    1158           7 :         if (!a) return 1;
    1159           7 :         b = NULL;
    1160             :       }
    1161             :       /* a or b is != NULL iff it is non-zero t_REAL; one of them is */
    1162         455 :       e = gexpo(y);
    1163         455 :       return (!a || r_approx0(a, e)) && (!b || r_approx0(b, e));
    1164     1748386 :     case t_REAL:
    1165     1748386 :       return !signe(x) || r_approx0(x, gexpo(y));
    1166      728674 :     default:
    1167      728674 :       return gequal0(x);
    1168             :   }
    1169             : }
    1170             : /*******************************************************************/
    1171             : /*                                                                 */
    1172             : /*                          VALUATION                              */
    1173             : /*             p is either a t_INT or a t_POL.                     */
    1174             : /*  returns the largest exponent of p dividing x when this makes   */
    1175             : /*  sense : error for types real, integermod and polymod if p does */
    1176             : /*  not divide the modulus, q-adic if q!=p.                        */
    1177             : /*                                                                 */
    1178             : /*******************************************************************/
    1179             : 
    1180             : static long
    1181         336 : minval(GEN x, GEN p)
    1182             : {
    1183         336 :   long i,k, val = LONG_MAX, lx = lg(x);
    1184        6538 :   for (i=lontyp[typ(x)]; i<lx; i++)
    1185             :   {
    1186        6202 :     k = gvaluation(gel(x,i),p);
    1187        6202 :     if (k < val) val = k;
    1188             :   }
    1189         336 :   return val;
    1190             : }
    1191             : 
    1192             : static int
    1193          91 : intdvd(GEN x, GEN y, GEN *z) { GEN r; *z = dvmdii(x,y,&r); return (r==gen_0); }
    1194             : 
    1195             : /* x t_FRAC, p t_INT, return v_p(x) */
    1196             : static long
    1197      292590 : frac_val(GEN x, GEN p) {
    1198      292590 :   long v = Z_pval(gel(x,2),p);
    1199      292590 :   if (v) return -v;
    1200      292454 :   return Z_pval(gel(x,1),p);
    1201             : }
    1202             : long
    1203     9336536 : Q_pval(GEN x, GEN p)
    1204             : {
    1205     9336536 :   if (lgefint(p) == 3) return Q_lval(x, uel(p,2));
    1206         552 :   return (typ(x)==t_INT)? Z_pval(x, p): frac_val(x, p);
    1207             : }
    1208             : 
    1209             : static long
    1210      392199 : frac_lval(GEN x, ulong p) {
    1211      392199 :   long v = Z_lval(gel(x,2),p);
    1212      392200 :   if (v) return -v;
    1213      222082 :   return Z_lval(gel(x,1),p);
    1214             : }
    1215             : long
    1216     9340731 : Q_lval(GEN x, ulong p){return (typ(x)==t_INT)? Z_lval(x, p): frac_lval(x, p);}
    1217             : 
    1218             : long
    1219     6304240 : Q_pvalrem(GEN x, GEN p, GEN *y)
    1220             : {
    1221             :   GEN a, b;
    1222             :   long v;
    1223     6304240 :   if (lgefint(p) == 3) return Q_lvalrem(x, uel(p,2), y);
    1224        5809 :   if (typ(x) == t_INT) return Z_pvalrem(x, p, y);
    1225           0 :   a = gel(x,1);
    1226           0 :   b = gel(x,2);
    1227           0 :   v = Z_pvalrem(b, p, &b);
    1228           0 :   if (v) { *y = isint1(b)? a: mkfrac(a, b); return -v; }
    1229           0 :   v = Z_pvalrem(a, p, &a);
    1230           0 :   *y = mkfrac(a, b); return v;
    1231             : }
    1232             : long
    1233     6302644 : Q_lvalrem(GEN x, ulong p, GEN *y)
    1234             : {
    1235             :   GEN a, b;
    1236             :   long v;
    1237     6302644 :   if (typ(x) == t_INT) return Z_lvalrem(x, p, y);
    1238      398340 :   a = gel(x,1);
    1239      398340 :   b = gel(x,2);
    1240      398340 :   v = Z_lvalrem(b, p, &b);
    1241      398343 :   if (v) { *y = isint1(b)? a: mkfrac(a, b); return -v; }
    1242      228917 :   v = Z_lvalrem(a, p, &a);
    1243      228917 :   *y = mkfrac(a, b); return v;
    1244             : }
    1245             : 
    1246             : long
    1247     1166401 : gvaluation(GEN x, GEN p)
    1248             : {
    1249     1166401 :   long tx = typ(x), tp;
    1250             :   pari_sp av;
    1251             : 
    1252     1166401 :   if (!p)
    1253          28 :     switch(tx)
    1254             :     {
    1255           7 :       case t_PADIC: return valp(x);
    1256           7 :       case t_POL: return RgX_val(x);
    1257           7 :       case t_SER: return valser(x);
    1258           7 :       default: pari_err_TYPE("gvaluation", x);
    1259             :     }
    1260     1166373 :   tp  = typ(p);
    1261     1166373 :   switch(tp)
    1262             :   {
    1263     1159359 :     case t_INT:
    1264     1159359 :       if (signe(p) && !is_pm1(p)) break;
    1265          28 :       pari_err_DOMAIN("gvaluation", "p", "=", p, p);
    1266        7007 :     case t_POL:
    1267        7007 :       if (degpol(p) > 0) break;
    1268             :     default:
    1269           7 :       pari_err_DOMAIN("gvaluation", "p", "=", p, p);
    1270             :   }
    1271             : 
    1272     1166338 :   switch(tx)
    1273             :   {
    1274      145061 :     case t_INT:
    1275      145061 :       if (!signe(x)) return LONG_MAX;
    1276      144956 :       if (tp == t_POL) return 0;
    1277      144634 :       return Z_pval(x,p);
    1278             : 
    1279          49 :     case t_REAL:
    1280          49 :       if (tp == t_POL) return 0;
    1281          21 :       break;
    1282             : 
    1283          28 :     case t_FFELT:
    1284          28 :       if (tp == t_POL) return FF_equal0(x)? LONG_MAX: 0;
    1285          14 :       break;
    1286             : 
    1287         105 :     case t_INTMOD: {
    1288         105 :       GEN a = gel(x,1), b = gel(x,2);
    1289             :       long val;
    1290         133 :       if (tp == t_POL) return signe(b)? 0: LONG_MAX;
    1291          42 :       av = avma;
    1292          42 :       if (!intdvd(a, p, &a)) break;
    1293          28 :       if (!intdvd(b, p, &b)) return gc_long(av,0);
    1294          14 :       val = 1; while (intdvd(a,p,&a) && intdvd(b,p,&b)) val++;
    1295          14 :       return gc_long(av,val);
    1296             :     }
    1297             : 
    1298      292509 :     case t_FRAC:
    1299      292509 :       if (tp == t_POL) return 0;
    1300      292495 :       return frac_val(x, p);
    1301             : 
    1302      721817 :     case t_PADIC:
    1303      721817 :       if (tp == t_POL) return 0;
    1304      721796 :       if (!equalii(p,gel(x,2))) break;
    1305      721789 :       return valp(x);
    1306             : 
    1307          35 :     case t_POLMOD: {
    1308          35 :       GEN a = gel(x,1), b = gel(x,2);
    1309             :       long v, val;
    1310          35 :       if (tp == t_INT) return gvaluation(b,p);
    1311          21 :       v = varn(p);
    1312          21 :       if (varn(a) != v) return 0;
    1313          21 :       av = avma;
    1314          21 :       a = RgX_divrem(a, p, ONLY_DIVIDES);
    1315          21 :       if (!a) break;
    1316          28 :       if (typ(b) != t_POL || varn(b) != v ||
    1317          21 :           !(b = RgX_divrem(b, p, ONLY_DIVIDES)) ) return gc_long(av,0);
    1318           7 :       val = 1;
    1319          28 :       while ((a = RgX_divrem(a, p, ONLY_DIVIDES)) &&
    1320          21 :              (b = RgX_divrem(b, p, ONLY_DIVIDES)) ) val++;
    1321           7 :       return gc_long(av,val);
    1322             :     }
    1323        6069 :     case t_POL: {
    1324        6069 :       if (tp == t_POL) {
    1325        5880 :         long vp = varn(p), vx = varn(x);
    1326        5880 :         if (vp == vx)
    1327             :         {
    1328             :           long val;
    1329        5866 :           if (RgX_is_monomial(p))
    1330             :           {
    1331        5831 :             val = RgX_val(x); if (val == LONG_MAX) return LONG_MAX;
    1332        5754 :             return val / degpol(p);
    1333             :           }
    1334          35 :           if (!signe(x)) return LONG_MAX;
    1335          21 :           av = avma;
    1336          21 :           for (val=0; ; val++)
    1337             :           {
    1338          35 :             x = RgX_divrem(x,p,ONLY_DIVIDES);
    1339          35 :             if (!x) return gc_long(av,val);
    1340          14 :             if (gc_needed(av,1))
    1341             :             {
    1342           0 :               if(DEBUGMEM>1) pari_warn(warnmem,"gvaluation");
    1343           0 :               x = gerepilecopy(av, x);
    1344             :             }
    1345             :           }
    1346             :         }
    1347          14 :         if (varncmp(vx, vp) > 0) return 0;
    1348             :       }
    1349         196 :       return minval(x,p);
    1350             :     }
    1351             : 
    1352         490 :     case t_SER: {
    1353         490 :       if (tp == t_POL) {
    1354         476 :         long vp = varn(p), vx = varn(x);
    1355         476 :         if (vp == vx)
    1356             :         {
    1357         469 :           long val = RgX_val(p);
    1358         469 :           if (!val) pari_err_DOMAIN("gvaluation", "p", "=", p, p);
    1359         462 :           return (long)(valser(x) / val);
    1360             :         }
    1361           7 :         if (varncmp(vx, vp) > 0) return 0;
    1362             :       }
    1363          14 :       return minval(x,p);
    1364             :     }
    1365             : 
    1366          49 :     case t_RFRAC:
    1367          49 :       return gvaluation(gel(x,1),p) - gvaluation(gel(x,2),p);
    1368             : 
    1369         126 :     case t_COMPLEX: case t_QUAD: case t_VEC: case t_COL: case t_MAT:
    1370         126 :       return minval(x,p);
    1371             :   }
    1372          63 :   pari_err_OP("valuation", x,p);
    1373             :   return 0; /* LCOV_EXCL_LINE */
    1374             : }
    1375             : GEN
    1376        3808 : gpvaluation(GEN x, GEN p)
    1377             : {
    1378        3808 :   long v = gvaluation(x,p);
    1379        3703 :   return v == LONG_MAX? mkoo(): stoi(v);
    1380             : }
    1381             : 
    1382             : /* x is nonzero */
    1383             : long
    1384    78023575 : u_lvalrem(ulong x, ulong p, ulong *py)
    1385             : {
    1386             :   ulong vx;
    1387    78023575 :   if (p == 2) { vx = vals(x); *py = x >> vx; return vx; }
    1388    68030996 :   for(vx = 0;;)
    1389             :   {
    1390   118967057 :     if (x % p) { *py = x; return vx; }
    1391    50936061 :     x /= p; /* gcc is smart enough to make a single div */
    1392    50936061 :     vx++;
    1393             :   }
    1394             : }
    1395             : long
    1396    65809901 : u_lval(ulong x, ulong p)
    1397             : {
    1398             :   ulong vx;
    1399    65809901 :   if (p == 2) return vals(x);
    1400    62523013 :   for(vx = 0;;)
    1401             :   {
    1402   101720294 :     if (x % p) return vx;
    1403    39197281 :     x /= p; /* gcc is smart enough to make a single div */
    1404    39197281 :     vx++;
    1405             :   }
    1406             : }
    1407             : 
    1408             : long
    1409     1824621 : z_lval(long s, ulong p) { return u_lval(labs(s), p); }
    1410             : long
    1411       87341 : z_lvalrem(long s, ulong p, long *py)
    1412             : {
    1413             :   long v;
    1414       87341 :   if (s < 0)
    1415             :   {
    1416           0 :     ulong u = (ulong)-s;
    1417           0 :     v = u_lvalrem(u, p, &u);
    1418           0 :     *py = -(long)u;
    1419             :   }
    1420             :   else
    1421             :   {
    1422       87341 :     ulong u = (ulong)s;
    1423       87341 :     v = u_lvalrem(u, p, &u);
    1424       87341 :     *py = (long)u;
    1425             :   }
    1426       87341 :   return v;
    1427             : }
    1428             : /* assume |p| > 1 */
    1429             : long
    1430     1318010 : z_pval(long s, GEN p)
    1431             : {
    1432     1318010 :   if (lgefint(p) > 3) return 0;
    1433     1318010 :   return z_lval(s, uel(p,2));
    1434             : }
    1435             : /* assume |p| > 1 */
    1436             : long
    1437         399 : z_pvalrem(long s, GEN p, long *py)
    1438             : {
    1439         399 :   if (lgefint(p) > 3) { *py = s; return 0; }
    1440         399 :   return z_lvalrem(s, uel(p,2), py);
    1441             : }
    1442             : 
    1443             : /* return v_q(x) and set *py = x / q^v_q(x), using divide & conquer */
    1444             : static long
    1445     2104787 : Z_pvalrem_DC(GEN x, GEN q, GEN *py)
    1446             : {
    1447     2104787 :   GEN r, z = dvmdii(x, q, &r);
    1448             :   long v;
    1449     2104781 :   if (r != gen_0) { *py = x; return 0; }
    1450     1452285 :   if (2 * lgefint(q) <= lgefint(z)+3) /* avoid squaring if pointless */
    1451     1436026 :     v = Z_pvalrem_DC(z, sqri(q), py) << 1;
    1452             :   else
    1453       16259 :   { v = 0; *py = z; }
    1454     1452286 :   z = dvmdii(*py, q, &r);
    1455     1452297 :   if (r != gen_0) return v + 1;
    1456      604751 :   *py = z; return v + 2;
    1457             : }
    1458             : 
    1459             : static const long VAL_DC_THRESHOLD = 16;
    1460             : 
    1461             : long
    1462    62778124 : Z_lval(GEN x, ulong p)
    1463             : {
    1464             :   long vx;
    1465             :   pari_sp av;
    1466    62778124 :   if (p == 2) return vali(x);
    1467    48410345 :   if (lgefint(x) == 3) return u_lval(uel(x,2), p);
    1468     2183952 :   av = avma;
    1469     2183952 :   for(vx = 0;;)
    1470    10589526 :   {
    1471             :     ulong r;
    1472    12773478 :     GEN q = absdiviu_rem(x, p, &r);
    1473    12774167 :     if (r) break;
    1474    10767327 :     vx++; x = q;
    1475    10767327 :     if (vx == VAL_DC_THRESHOLD) {
    1476      177801 :       if (p == 1) pari_err_DOMAIN("Z_lval", "p", "=", gen_1, gen_1);
    1477      177801 :       vx += Z_pvalrem_DC(x, sqru(p), &x) << 1;
    1478      177801 :       q = absdiviu_rem(x, p, &r); if (!r) vx++;
    1479      177801 :       break;
    1480             :     }
    1481             :   }
    1482     2184641 :   return gc_long(av,vx);
    1483             : }
    1484             : long
    1485    59611957 : Z_lvalrem(GEN x, ulong p, GEN *py)
    1486             : {
    1487             :   long vx, sx;
    1488             :   pari_sp av;
    1489    59611957 :   if (p == 2) { vx = vali(x); *py = shifti(x, -vx); return vx; }
    1490    46171126 :   if (lgefint(x) == 3) {
    1491             :     ulong u;
    1492    40064558 :     vx = u_lvalrem(uel(x,2), p, &u);
    1493    40064216 :     *py = signe(x) < 0? utoineg(u): utoipos(u);
    1494    40063077 :     return vx;
    1495             :   }
    1496     6106568 :   av = avma; (void)new_chunk(lgefint(x));
    1497     6106870 :   sx = signe(x);
    1498     6106870 :   for(vx = 0;;)
    1499    16521980 :   {
    1500             :     ulong r;
    1501    22628850 :     GEN q = absdiviu_rem(x, p, &r);
    1502    22628842 :     if (r) break;
    1503    17012345 :     vx++; x = q;
    1504    17012345 :     if (vx == VAL_DC_THRESHOLD) {
    1505      490365 :       if (p == 1) pari_err_DOMAIN("Z_lvalrem", "p", "=", gen_1, gen_1);
    1506      490365 :       vx += Z_pvalrem_DC(x, sqru(p), &x) << 1;
    1507      490363 :       q = absdiviu_rem(x, p, &r); if (!r) { vx++; x = q; }
    1508      490362 :       break;
    1509             :     }
    1510             :   }
    1511     6106859 :   set_avma(av); *py = icopy(x); setsigne(*py, sx); return vx;
    1512             : }
    1513             : 
    1514             : /* Is |q| <= p ? */
    1515             : static int
    1516    21938380 : isless_iu(GEN q, ulong p) {
    1517    21938380 :   long l = lgefint(q);
    1518    21938380 :   return l==2 || (l == 3 && uel(q,2) <= p);
    1519             : }
    1520             : 
    1521             : long
    1522   247379456 : u_lvalrem_stop(ulong *n, ulong p, int *stop)
    1523             : {
    1524   247379456 :   ulong N = *n, q = N / p, r = N % p; /* gcc makes a single div */
    1525   247379456 :   long v = 0;
    1526   247379456 :   if (!r)
    1527             :   {
    1528    29248093 :     do { v++; N = q; q = N / p; r = N % p; } while (!r);
    1529    20572969 :     *n = N;
    1530             :   }
    1531   247379456 :   *stop = q <= p; return v;
    1532             : }
    1533             : /* Assume n > 0. Return v_p(n), set *n := n/p^v_p(n). Set 'stop' if now
    1534             :  * n < p^2 [implies n prime if no prime < p divides n] */
    1535             : long
    1536   155890494 : Z_lvalrem_stop(GEN *n, ulong p, int *stop)
    1537             : {
    1538             :   pari_sp av;
    1539             :   long v;
    1540             :   ulong r;
    1541             :   GEN N, q;
    1542             : 
    1543   155890494 :   if (lgefint(*n) == 3)
    1544             :   {
    1545   133952856 :     r = (*n)[2];
    1546   133952856 :     v = u_lvalrem_stop(&r, p, stop);
    1547   133954021 :     if (v) *n = utoipos(r);
    1548   133952355 :     return v;
    1549             :   }
    1550    21937638 :   av = avma; v = 0; q = absdiviu_rem(*n, p, &r);
    1551    21938392 :   if (r) set_avma(av);
    1552             :   else
    1553             :   {
    1554             :     do {
    1555      267210 :       v++; N = q;
    1556      267210 :       if (v == VAL_DC_THRESHOLD)
    1557             :       {
    1558         607 :         v += Z_pvalrem_DC(N,sqru(p),&N) << 1;
    1559         607 :         q = absdiviu_rem(N, p, &r); if (!r) { v++; N = q; }
    1560         607 :         break;
    1561             :       }
    1562      266603 :       q = absdiviu_rem(N, p, &r);
    1563      266603 :     } while (!r);
    1564      226158 :     *n = N;
    1565             :   }
    1566    21938395 :   *stop = isless_iu(q,p); return v;
    1567             : }
    1568             : 
    1569             : /* x is a nonzero integer, |p| > 1 */
    1570             : long
    1571    49685245 : Z_pvalrem(GEN x, GEN p, GEN *py)
    1572             : {
    1573             :   long vx;
    1574             :   pari_sp av;
    1575             : 
    1576    49685245 :   if (lgefint(p) == 3) return Z_lvalrem(x, uel(p,2), py);
    1577       25689 :   if (lgefint(x) == 3) { *py = icopy(x); return 0; }
    1578        8355 :   av = avma; vx = 0; (void)new_chunk(lgefint(x));
    1579             :   for(;;)
    1580       14251 :   {
    1581       22672 :     GEN r, q = dvmdii(x,p,&r);
    1582       22672 :     if (r != gen_0) { set_avma(av); *py = icopy(x); return vx; }
    1583       14251 :     vx++; x = q;
    1584             :   }
    1585             : }
    1586             : long
    1587     2595626 : u_pvalrem(ulong x, GEN p, ulong *py)
    1588             : {
    1589     2595626 :   if (lgefint(p) == 3) return u_lvalrem(x, uel(p,2), py);
    1590         549 :   *py = x; return 0;
    1591             : }
    1592             : long
    1593      136960 : u_pval(ulong x, GEN p)
    1594             : {
    1595      136960 :   if (lgefint(p) == 3) return u_lval(x, uel(p,2));
    1596           0 :   return 0;
    1597             : }
    1598             : long
    1599    47061461 : Z_pval(GEN x, GEN p) {
    1600             :   long vx;
    1601             :   pari_sp av;
    1602             : 
    1603    47061461 :   if (lgefint(p) == 3) return Z_lval(x, uel(p,2));
    1604       32938 :   if (lgefint(x) == 3) return 0;
    1605        7961 :   av = avma; vx = 0;
    1606             :   for(;;)
    1607       25450 :   {
    1608       33411 :     GEN r, q = dvmdii(x,p,&r);
    1609       33450 :     if (r != gen_0) return gc_long(av,vx);
    1610       25450 :     vx++; x = q;
    1611             :   }
    1612             : }
    1613             : 
    1614             : /* return v_p(n!) = [n/p] + [n/p^2] + ... */
    1615             : long
    1616     1986281 : factorial_lval(ulong n, ulong p)
    1617             : {
    1618             :   ulong q, v;
    1619     1986281 :   if (p == 2) return n - hammingl(n);
    1620     1317732 :   q = p; v = 0;
    1621     1446274 :   do { v += n/q; q *= p; } while (n >= q);
    1622     1317732 :   return (long)v;
    1623             : }
    1624             : 
    1625             : /********** Same for "containers" ZX / ZV / ZC **********/
    1626             : 
    1627             : /* If the t_INT q divides the ZX/ZV x, return the quotient. Otherwise NULL.
    1628             :  * Stack clean; assumes lg(x) > 1 */
    1629             : static GEN
    1630        6580 : gen_Z_divides(GEN x, GEN q, long imin)
    1631             : {
    1632             :   long i, l;
    1633        6580 :   GEN y = cgetg_copy(x, &l);
    1634             : 
    1635        6580 :   y[1] = x[1]; /* Needed for ZX; no-op if ZV, overwritten in first iteration */
    1636       88544 :   for (i = imin; i < l; i++)
    1637             :   {
    1638       85960 :     GEN r, xi = gel(x,i);
    1639       85960 :     if (!signe(xi)) { gel(y,i) = xi; continue; }
    1640       55390 :     gel(y,i) = dvmdii(xi, q, &r);
    1641       55390 :     if (r != gen_0) { set_avma((pari_sp)(y+l)); return NULL; }
    1642             :   }
    1643        2584 :   return y;
    1644             : }
    1645             : /* If q divides the ZX/ZV x, return the quotient. Otherwise NULL.
    1646             :  * Stack clean; assumes lg(x) > 1 */
    1647             : static GEN
    1648        4885 : gen_z_divides(GEN x, ulong q, long imin)
    1649             : {
    1650             :   long i, l;
    1651        4885 :   GEN y = cgetg_copy(x, &l);
    1652             : 
    1653        4885 :   y[1] = x[1]; /* Needed for ZX; no-op if ZV, overwritten in first iteration */
    1654       43230 :   for (i = imin; i < l; i++)
    1655             :   {
    1656             :     ulong r;
    1657       41929 :     GEN xi = gel(x,i);
    1658       41929 :     if (!signe(xi)) { gel(y,i) = xi; continue; }
    1659       28615 :     gel(y,i) = absdiviu_rem(xi, q, &r);
    1660       28615 :     if (r) { set_avma((pari_sp)(y+l)); return NULL; }
    1661       25031 :     affectsign_safe(xi, &gel(y,i));
    1662             :   }
    1663        1301 :   return y;
    1664             : }
    1665             : 
    1666             : /* return v_q(x) and set *py = x / q^v_q(x), using divide & conquer */
    1667             : static long
    1668       11428 : gen_pvalrem_DC(GEN x, GEN q, GEN *py, long imin)
    1669             : {
    1670             : 
    1671       11428 :   pari_sp av = avma;
    1672       11428 :   long v, i, l, lz = LONG_MAX;
    1673       11428 :   GEN y = cgetg_copy(x, &l);
    1674             : 
    1675       11428 :   y[1] = x[1];
    1676      134570 :   for (i = imin; i < l; i++)
    1677             :   {
    1678      127990 :     GEN r, xi = gel(x,i);
    1679      127990 :     if (!signe(xi)) { gel(y,i) = xi; continue; }
    1680       86178 :     gel(y,i) = dvmdii(xi, q, &r);
    1681       86178 :     if (r != gen_0) { *py = x; return gc_long(av,0); }
    1682       81330 :     lz = minss(lz, lgefint(gel(y,i)));
    1683             :   }
    1684        6580 :   if (2 * lgefint(q) <= lz+3) /* avoid squaring if pointless */
    1685        6529 :     v = gen_pvalrem_DC(y, sqri(q), py, imin) << 1;
    1686             :   else
    1687          51 :   { v = 0; *py = y; }
    1688             : 
    1689        6580 :   y = gen_Z_divides(*py, q, imin);
    1690        6580 :   if (!y) return v+1;
    1691        2584 :   *py = y; return v+2;
    1692             : }
    1693             : 
    1694             : static long
    1695      775319 : gen_2val(GEN x, long imin)
    1696             : {
    1697      775319 :   long i, lx = lg(x), v = LONG_MAX;
    1698     2925980 :   for (i = imin; i < lx; i++)
    1699             :   {
    1700     2488536 :     GEN c = gel(x,i);
    1701             :     long w;
    1702     2488536 :     if (!signe(c)) continue;
    1703     2275615 :     w = vali(c);
    1704     2275615 :     if (w < v) { v = w; if (!v) break; }
    1705             :   }
    1706      775319 :   return v;
    1707             : }
    1708             : static long
    1709      869811 : gen_lval(GEN x, ulong p, long imin)
    1710             : {
    1711             :   long i, lx, v;
    1712             :   pari_sp av;
    1713             :   GEN y;
    1714      869811 :   if (p == 2) return gen_2val(x, imin);
    1715       94492 :   av = avma;
    1716       94492 :   lx = lg(x); y = leafcopy(x);
    1717      274782 :   for(v = 0;; v++)
    1718     1376795 :     for (i = imin; i < lx; i++)
    1719             :     {
    1720             :       ulong r;
    1721     1196505 :       gel(y,i) = absdiviu_rem(gel(y,i), p, &r);
    1722     1196505 :       if (r) return gc_long(av,v);
    1723             :     }
    1724             : }
    1725             : long
    1726      748784 : ZX_lval(GEN x, ulong p) { return gen_lval(x, p, 2); }
    1727             : long
    1728           0 : ZV_lval(GEN x, ulong p) { return gen_lval(x, p, 1); }
    1729             : 
    1730             : long
    1731       28931 : zx_lval(GEN f, long p)
    1732             : {
    1733       28931 :   long i, l = lg(f), x = LONG_MAX;
    1734       30219 :   for(i=2; i<l; i++)
    1735             :   {
    1736             :     long y;
    1737       29477 :     if (f[i] == 0) continue;
    1738       29428 :     y = z_lval(f[i], p);
    1739       29428 :     if (y < x) { x = y; if (x == 0) return x; }
    1740             :   }
    1741         742 :   return x;
    1742             : }
    1743             : 
    1744             : static long
    1745      131345 : gen_pval(GEN x, GEN p, long imin)
    1746             : {
    1747             :   long i, lx, v;
    1748             :   pari_sp av;
    1749             :   GEN y;
    1750      131345 :   if (lgefint(p) == 3) return gen_lval(x, p[2], imin);
    1751       10318 :   av = avma;
    1752       10318 :   lx = lg(x); y = leafcopy(x);
    1753       10318 :   for(v = 0;; v++)
    1754             :   {
    1755       10318 :     if (v == VAL_DC_THRESHOLD)
    1756             :     {
    1757           0 :       if (is_pm1(p)) pari_err_DOMAIN("gen_pval", "p", "=", p, p);
    1758           0 :       v += gen_pvalrem_DC(y, p, &y, imin);
    1759           0 :       return gc_long(av,v);
    1760             :     }
    1761             : 
    1762       10318 :     for (i = imin; i < lx; i++)
    1763             :     {
    1764       10318 :       GEN r; gel(y,i) = dvmdii(gel(y,i), p, &r);
    1765       10318 :       if (r != gen_0) return gc_long(av,v);
    1766             :     }
    1767             :   }
    1768             : }
    1769             : long
    1770       97705 : ZX_pval(GEN x, GEN p) { return gen_pval(x, p, 2); }
    1771             : long
    1772       33640 : ZV_pval(GEN x, GEN p) { return gen_pval(x, p, 1); }
    1773             : /* v = 0 (mod p) */
    1774             : int
    1775         616 : ZV_Z_dvd(GEN v, GEN p)
    1776             : {
    1777         616 :   pari_sp av = avma;
    1778         616 :   long i, l = lg(v);
    1779        1911 :   for (i=1; i<l; i++)
    1780        1407 :     if (!dvdii(gel(v,i), p)) return gc_long(av,0);
    1781         504 :   return gc_long(av,1);
    1782             : }
    1783             : 
    1784             : static long
    1785     4801995 : gen_2valrem(GEN x, GEN *px, long imin)
    1786             : {
    1787     4801995 :   long i, lx = lg(x), v = LONG_MAX;
    1788             :   GEN z;
    1789    13777387 :   for (i = imin; i < lx; i++)
    1790             :   {
    1791    12465828 :     GEN c = gel(x,i);
    1792             :     long w;
    1793    12465828 :     if (!signe(c)) continue;
    1794    11582080 :     w = vali(c);
    1795    11582307 :     if (w < v) {
    1796     6897082 :       v = w;
    1797     6897082 :       if (!v) { *px = x; return 0; } /* early abort */
    1798             :     }
    1799             :   }
    1800     1311559 :   z = cgetg_copy(x, &lx); z[1] = x[1];
    1801     8412063 :   for (i=imin; i<lx; i++) gel(z,i) = shifti(gel(x,i), -v);
    1802     1311290 :   *px = z; return v;
    1803             : }
    1804             : static long
    1805     8209355 : gen_lvalrem(GEN x, ulong p, GEN *px, long imin)
    1806             : {
    1807             :   long i, lx, v;
    1808             :   GEN y;
    1809     8209355 :   if (p == 2) return gen_2valrem(x, px, imin);
    1810     3407414 :   y = cgetg_copy(x, &lx);
    1811     3407678 :   y[1] = x[1];
    1812     3407678 :   x = leafcopy(x);
    1813     3407156 :   for(v = 0;; v++)
    1814             :   {
    1815     4789651 :     if (v == VAL_DC_THRESHOLD)
    1816             :     {
    1817        4885 :       if (p == 1) pari_err_DOMAIN("gen_lvalrem", "p", "=", gen_1, gen_1);
    1818        4885 :       v += gen_pvalrem_DC(x, sqru(p), px, imin) << 1;
    1819        4885 :       x = gen_z_divides(*px, p, imin);
    1820        4885 :       if (x) { *px = x; v++; }
    1821        4885 :       return v;
    1822             :     }
    1823             : 
    1824    15277054 :     for (i = imin; i < lx; i++)
    1825             :     {
    1826    13894559 :       ulong r; gel(y,i) = absdiviu_rem(gel(x,i), p, &r);
    1827    13894348 :       if (r) { *px = x; return v; }
    1828    10491609 :       affectsign_safe(gel(x,i), &gel(y,i));
    1829             :     }
    1830     1382495 :     swap(x, y);
    1831             :   }
    1832             : }
    1833             : long
    1834         721 : ZX_lvalrem(GEN x, ulong p, GEN *px) { return gen_lvalrem(x,p,px, 2); }
    1835             : long
    1836           0 : ZV_lvalrem(GEN x, ulong p, GEN *px) { return gen_lvalrem(x,p,px, 1); }
    1837             : 
    1838             : static long
    1839     8221543 : gen_pvalrem(GEN x, GEN p, GEN *px, long imin)
    1840             : {
    1841             :   long i, lx, v;
    1842             :   GEN y;
    1843     8221543 :   if (lgefint(p) == 3) return gen_lvalrem(x, p[2], px, imin);
    1844       12943 :   y = cgetg_copy(x, &lx);
    1845       12945 :   y[1] = x[1];
    1846       12945 :   x = leafcopy(x);
    1847       12945 :   for(v = 0;; v++)
    1848             :   {
    1849       13734 :     if (v == VAL_DC_THRESHOLD)
    1850             :     {
    1851          14 :       if (is_pm1(p)) pari_err_DOMAIN("gen_pvalrem", "p", "=", p, p);
    1852          14 :       return v + gen_pvalrem_DC(x, p, px, imin);
    1853             :     }
    1854             : 
    1855       22384 :     for (i = imin; i < lx; i++)
    1856             :     {
    1857       21595 :       GEN r; gel(y,i) = dvmdii(gel(x,i), p, &r);
    1858       21595 :       if (r != gen_0) { *px = x; return v; }
    1859             :     }
    1860         789 :     swap(x, y);
    1861             :   }
    1862             : }
    1863             : long
    1864     4252779 : ZX_pvalrem(GEN x, GEN p, GEN *px) { return gen_pvalrem(x,p,px, 2); }
    1865             : long
    1866     3968766 : ZV_pvalrem(GEN x, GEN p, GEN *px) { return gen_pvalrem(x,p,px, 1); }
    1867             : 
    1868             : static long
    1869           0 : ZX_gen_pvalrem(GEN x, GEN p, GEN *px, long imin)
    1870             : {
    1871             :   long i, lx, v;
    1872             :   GEN y;
    1873           0 :   y = cgetg_copy(x, &lx);
    1874           0 :   y[1] = x[1];
    1875           0 :   x = leafcopy(x);
    1876           0 :   for (i = imin; i < lx; i++)
    1877           0 :     if (typ(gel(x, i)) != t_INT)
    1878             :     {
    1879           0 :       gel(x, i) = leafcopy(gel(x,i));
    1880           0 :       gel(y, i) = leafcopy(gel(x,i));
    1881             :     }
    1882           0 :   for(v = 0;; v++)
    1883             :   {
    1884             : #if 0
    1885             :     if (v == VAL_DC_THRESHOLD) /* TODO */
    1886             :     {
    1887             :       if (is_pm1(p)) pari_err_DOMAIN("ZX_gen_pvalrem", "p", "=", p, p);
    1888             :       return v + ZX_gen_pvalrem_DC(x, p, px, imin);
    1889             :     }
    1890             : #endif
    1891             : 
    1892           0 :     for (i = imin; i < lx; i++)
    1893             :     {
    1894           0 :       GEN r, xi = gel(x,i);
    1895           0 :       if (typ(xi) == t_INT)
    1896             :       {
    1897           0 :         gel(y,i) = dvmdii(xi, p, &r);
    1898           0 :         if (r != gen_0) { *px = x; return v; }
    1899             :       } else
    1900             :       {
    1901           0 :         long j, lxi = lg(xi);
    1902           0 :         for(j = 2; j < lxi; j++)
    1903             :         {
    1904           0 :           gmael(y,i,j) = dvmdii(gel(xi,j), p, &r);
    1905           0 :           if (r != gen_0) { *px = x; return v; }
    1906             :         }
    1907             :       }
    1908             :     }
    1909           0 :     swap(x, y);
    1910             :   }
    1911             : }
    1912             : 
    1913             : long
    1914           0 : ZXX_pvalrem(GEN x, GEN p, GEN *px) { return ZX_gen_pvalrem(x,p,px, 2); }
    1915             : long
    1916           0 : ZXV_pvalrem(GEN x, GEN p, GEN *px) { return ZX_gen_pvalrem(x,p,px, 1); }
    1917             : 
    1918             : /*******************************************************************/
    1919             : /*                                                                 */
    1920             : /*                       NEGATION: Create -x                       */
    1921             : /*                                                                 */
    1922             : /*******************************************************************/
    1923             : 
    1924             : GEN
    1925   466446623 : gneg(GEN x)
    1926             : {
    1927             :   GEN y;
    1928   466446623 :   switch(typ(x))
    1929             :   {
    1930   137789724 :     case t_INT:
    1931   137789724 :       return signe(x)? negi(x): gen_0;
    1932   241640760 :     case t_REAL:
    1933   241640760 :       return mpneg(x);
    1934             : 
    1935      222171 :     case t_INTMOD: y=cgetg(3,t_INTMOD);
    1936      222171 :       gel(y,1) = icopy(gel(x,1));
    1937      222171 :       gel(y,2) = signe(gel(x,2))? subii(gel(y,1),gel(x,2)): gen_0;
    1938      222170 :       break;
    1939             : 
    1940     3140175 :     case t_FRAC:
    1941     3140175 :       y = cgetg(3, t_FRAC);
    1942     3140175 :       gel(y,1) = negi(gel(x,1));
    1943     3140173 :       gel(y,2) = icopy(gel(x,2)); break;
    1944             : 
    1945    75682024 :     case t_COMPLEX:
    1946    75682024 :       y=cgetg(3, t_COMPLEX);
    1947    75683606 :       gel(y,1) = gneg(gel(x,1));
    1948    75685442 :       gel(y,2) = gneg(gel(x,2));
    1949    75685166 :       break;
    1950             : 
    1951      246190 :     case t_POLMOD:
    1952      246190 :       retmkpolmod(gneg(gel(x,2)), RgX_copy(gel(x,1)));
    1953             : 
    1954      153370 :     case t_RFRAC:
    1955      153370 :       y = cgetg(3, t_RFRAC);
    1956      153370 :       gel(y,1) = gneg(gel(x,1));
    1957      153370 :       gel(y,2) = RgX_copy(gel(x,2)); break;
    1958             : 
    1959      650294 :     case t_PADIC:
    1960      650294 :       if (!signe(gel(x,4))) return gcopy(x);
    1961      646675 :       y = cgetg(5, t_PADIC);
    1962      646675 :       y[1] = x[1];
    1963      646675 :       gel(y,2) = icopy(gel(x,2));
    1964      646674 :       gel(y,3) = icopy(gel(x,3));
    1965      646674 :       gel(y,4) = subii(gel(x,3),gel(x,4));
    1966      646673 :       break;
    1967             : 
    1968         133 :     case t_QUAD:
    1969         133 :       y=cgetg(4,t_QUAD);
    1970         133 :       gel(y,1) = ZX_copy(gel(x,1));
    1971         133 :       gel(y,2) = gneg(gel(x,2));
    1972         133 :       gel(y,3) = gneg(gel(x,3)); break;
    1973             : 
    1974       80424 :     case t_FFELT: return FF_neg(x);
    1975     6490656 :     case t_POL: return RgX_neg(x);
    1976       16688 :     case t_SER: pari_APPLY_ser_normalized(gneg(gel(x,i)));
    1977        1533 :     case t_VEC: return RgV_neg(x);
    1978      462964 :     case t_COL: return RgC_neg(x);
    1979         315 :     case t_MAT: return RgM_neg(x);
    1980         777 :     case t_INFINITY: return inf_get_sign(x) == 1? mkmoo(): mkoo();
    1981           0 :     default:
    1982           0 :       pari_err_TYPE("gneg",x);
    1983             :       return NULL; /* LCOV_EXCL_LINE */
    1984             :   }
    1985    79847639 :   return y;
    1986             : }
    1987             : 
    1988             : GEN
    1989   133310537 : gneg_i(GEN x)
    1990             : {
    1991             :   GEN y;
    1992   133310537 :   switch(typ(x))
    1993             :   {
    1994    67789965 :     case t_INT:
    1995    67789965 :       return signe(x)? negi(x): gen_0;
    1996    31815626 :     case t_REAL:
    1997    31815626 :       return mpneg(x);
    1998             : 
    1999      717747 :     case t_INTMOD: y=cgetg(3,t_INTMOD);
    2000      717747 :       gel(y,1) = gel(x,1);
    2001      717747 :       gel(y,2) = signe(gel(x,2))? subii(gel(y,1),gel(x,2)): gen_0;
    2002      717747 :       break;
    2003             : 
    2004     5419392 :     case t_FRAC:
    2005     5419392 :       y = cgetg(3, t_FRAC);
    2006     5419389 :       gel(y,1) = negi(gel(x,1));
    2007     5419386 :       gel(y,2) = gel(x,2); break;
    2008             : 
    2009    10751570 :     case t_COMPLEX:
    2010    10751570 :       y = cgetg(3, t_COMPLEX);
    2011    10751690 :       gel(y,1) = gneg_i(gel(x,1));
    2012    10751737 :       gel(y,2) = gneg_i(gel(x,2)); break;
    2013             : 
    2014     2015113 :     case t_PADIC: y = cgetg(5,t_PADIC);
    2015     2015113 :       y[1] = x[1];
    2016     2015113 :       gel(y,2) = gel(x,2);
    2017     2015113 :       gel(y,3) = gel(x,3);
    2018     2015113 :       gel(y,4) = signe(gel(x,4))? subii(gel(x,3),gel(x,4)): gen_0; break;
    2019             : 
    2020      138177 :     case t_POLMOD:
    2021      138177 :       retmkpolmod(gneg_i(gel(x,2)), RgX_copy(gel(x,1)));
    2022             : 
    2023       84504 :     case t_FFELT: return FF_neg_i(x);
    2024             : 
    2025         672 :     case t_QUAD: y=cgetg(4,t_QUAD);
    2026         672 :       gel(y,1) = gel(x,1);
    2027         672 :       gel(y,2) = gneg_i(gel(x,2));
    2028         672 :       gel(y,3) = gneg_i(gel(x,3)); break;
    2029             : 
    2030        2478 :     case t_VEC:
    2031             :     case t_COL:
    2032       13594 :     case t_MAT: pari_APPLY_same(gneg_i(gel(x,i)));
    2033    36685402 :     case t_POL: pari_APPLY_pol_normalized(gneg_i(gel(x,i)));
    2034     3749049 :     case t_SER: pari_APPLY_ser_normalized(gneg_i(gel(x,i)));
    2035             : 
    2036     5044871 :     case t_RFRAC:
    2037     5044871 :       y = cgetg(3, t_RFRAC);
    2038     5044871 :       gel(y,1) = gneg_i(gel(x,1));
    2039     5044862 :       gel(y,2) = gel(x,2); break;
    2040             : 
    2041           0 :     default:
    2042           0 :       pari_err_TYPE("gneg_i",x);
    2043             :       return NULL; /* LCOV_EXCL_LINE */
    2044             :   }
    2045    23949507 :   return y;
    2046             : }
    2047             : 
    2048             : /******************************************************************/
    2049             : /*                                                                */
    2050             : /*                       ABSOLUTE VALUE                           */
    2051             : /*    Create abs(x) if x is integer, real, fraction or complex.   */
    2052             : /*                       Error otherwise.                         */
    2053             : /*                                                                */
    2054             : /******************************************************************/
    2055             : static int
    2056           0 : is_negative(GEN x) {
    2057           0 :   switch(typ(x))
    2058             :   {
    2059           0 :     case t_INT: case t_REAL:
    2060           0 :       return (signe(x) < 0);
    2061           0 :     case t_FRAC:
    2062           0 :       return (signe(gel(x,1)) < 0);
    2063             :   }
    2064           0 :   return 0;
    2065             : }
    2066             : 
    2067             : GEN
    2068    53203992 : gabs(GEN x, long prec)
    2069             : {
    2070             :   long lx;
    2071             :   pari_sp av;
    2072             :   GEN y, N;
    2073             : 
    2074    53203992 :   switch(typ(x))
    2075             :   {
    2076    35151844 :     case t_INT: case t_REAL:
    2077    35151844 :       return mpabs(x);
    2078             : 
    2079       12776 :     case t_FRAC:
    2080       12776 :       return absfrac(x);
    2081             : 
    2082    17943148 :     case t_COMPLEX:
    2083    17943148 :       av=avma; N=cxnorm(x);
    2084    17922868 :       switch(typ(N))
    2085             :       {
    2086         266 :         case t_INT:
    2087         266 :           if (!Z_issquareall(N, &y)) break;
    2088         105 :           return gerepileupto(av, y);
    2089       21735 :         case t_FRAC: {
    2090             :           GEN a,b;
    2091       36036 :           if (!Z_issquareall(gel(N,1), &a)) break;
    2092       14301 :           if (!Z_issquareall(gel(N,2), &b)) break;
    2093           0 :           return gerepileupto(av, gdiv(a,b));
    2094             :         }
    2095             :       }
    2096    17922763 :       return gerepileupto(av, gsqrt(N,prec));
    2097             : 
    2098          21 :     case t_QUAD:
    2099          21 :       av = avma;
    2100          21 :       return gerepileuptoleaf(av, gabs(quadtofp(x, prec), prec));
    2101             : 
    2102           0 :     case t_POL:
    2103           0 :       lx = lg(x); if (lx<=2) return RgX_copy(x);
    2104           0 :       return is_negative(gel(x,lx-1))? RgX_neg(x): RgX_copy(x);
    2105             : 
    2106           7 :     case t_SER:
    2107           7 :      if (!signe(x)) pari_err_DOMAIN("abs", "argument", "=", gen_0, x);
    2108           7 :      if (valser(x)) pari_err_DOMAIN("abs", "series valuation", "!=", gen_0, x);
    2109           0 :      return is_negative(gel(x,2))? gneg(x): gcopy(x);
    2110             : 
    2111      101963 :     case t_VEC: case t_COL: case t_MAT:
    2112      591394 :       pari_APPLY_same(gabs(gel(x,i),prec));
    2113             : 
    2114          14 :     case t_INFINITY:
    2115          14 :       return mkoo();
    2116             :   }
    2117           0 :   pari_err_TYPE("gabs",x);
    2118             :   return NULL; /* LCOV_EXCL_LINE */
    2119             : }
    2120             : 
    2121             : GEN
    2122       79310 : gmax(GEN x, GEN y) { return gcopy(gmax_shallow(x,y)); }
    2123             : GEN
    2124           0 : gmaxgs(GEN x, long s) { return (gcmpsg(s,x)>=0)? stoi(s): gcopy(x); }
    2125             : 
    2126             : GEN
    2127       12264 : gmin(GEN x, GEN y) { return gcopy(gmin_shallow(x,y)); }
    2128             : GEN
    2129           0 : gmings(GEN x, long s) { return (gcmpsg(s,x)>0)? gcopy(x): stoi(s); }
    2130             : 
    2131             : long
    2132      503192 : vecindexmax(GEN x)
    2133             : {
    2134      503192 :   long lx = lg(x), i0, i;
    2135             :   GEN s;
    2136             : 
    2137      503192 :   if (lx==1) pari_err_DOMAIN("vecindexmax", "empty argument", "=", x,x);
    2138      503194 :   switch(typ(x))
    2139             :   {
    2140      503194 :     case t_VEC: case t_COL:
    2141      503194 :       s = gel(x,i0=1);
    2142     1502464 :       for (i=2; i<lx; i++)
    2143      999268 :         if (gcmp(gel(x,i),s) > 0) s = gel(x,i0=i);
    2144      503196 :       return i0;
    2145           0 :     case t_VECSMALL:
    2146           0 :       return vecsmall_indexmax(x);
    2147           0 :     default: pari_err_TYPE("vecindexmax",x);
    2148             :   }
    2149             :   /* LCOV_EXCL_LINE */
    2150           0 :   return 0;
    2151             : }
    2152             : long
    2153      181366 : vecindexmin(GEN x)
    2154             : {
    2155      181366 :   long lx = lg(x), i0, i;
    2156             :   GEN s;
    2157             : 
    2158      181366 :   if (lx==1) pari_err_DOMAIN("vecindexmin", "empty argument", "=", x,x);
    2159      181366 :   switch(typ(x))
    2160             :   {
    2161      181366 :     case t_VEC: case t_COL:
    2162      181366 :       s = gel(x,i0=1);
    2163      943837 :       for (i=2; i<lx; i++)
    2164      762471 :         if (gcmp(gel(x,i),s) < 0) s = gel(x,i0=i);
    2165      181366 :       return i0;
    2166           0 :     case t_VECSMALL:
    2167           0 :       return vecsmall_indexmin(x);
    2168           0 :     default: pari_err_TYPE("vecindexmin",x);
    2169             :   }
    2170             :   /* LCOV_EXCL_LINE */
    2171           0 :   return 0;
    2172             : }
    2173             : 
    2174             : GEN
    2175      226753 : vecmax0(GEN x, GEN *pi)
    2176             : {
    2177      226753 :   long i, lx = lg(x), tx = typ(x);
    2178      226753 :   if (!is_matvec_t(tx) && tx != t_VECSMALL
    2179          49 :       && (tx != t_LIST || list_typ(x) != t_LIST_RAW)) return gcopy(x);
    2180      226732 :   if (tx == t_LIST)
    2181          28 :   { if (list_data(x)) { x = list_data(x); lx = lg(x); } else lx = 1; }
    2182      226732 :   if (lx==1) pari_err_DOMAIN("vecmax", "empty argument", "=", x,x);
    2183      226690 :   switch(typ(x))
    2184             :   {
    2185      226213 :     case t_VEC: case t_COL:
    2186      226213 :       i = vecindexmax(x); if (pi) *pi = utoipos(i);
    2187      226216 :       return gcopy(gel(x,i));
    2188         456 :     case t_MAT: {
    2189         456 :       long j, i0 = 1, j0 = 1, lx2 = lgcols(x);
    2190             :       GEN s;
    2191         456 :       if (lx2 == 1) pari_err_DOMAIN("vecmax", "empty argument", "=", x,x);
    2192         449 :       s = gcoeff(x,i0,j0); i = 2;
    2193        1205 :       for (j=1; j<lx; j++,i=1)
    2194             :       {
    2195         756 :         GEN c = gel(x,j);
    2196        1817 :         for (; i<lx2; i++)
    2197        1061 :           if (gcmp(gel(c,i),s) > 0) { s = gel(c,i); j0=j; i0=i; }
    2198             :       }
    2199         449 :       if (pi) *pi = mkvec2(utoipos(i0), utoipos(j0));
    2200         449 :       return gcopy(s);
    2201             :     }
    2202          21 :     case t_VECSMALL:
    2203          21 :       i = vecsmall_indexmax(x); if (pi) *pi = utoipos(i);
    2204          21 :       return stoi(x[i]);
    2205             :   }
    2206             :   return NULL;/*LCOV_EXCL_LINE*/
    2207             : }
    2208             : GEN
    2209      146716 : vecmin0(GEN x, GEN *pi)
    2210             : {
    2211      146716 :   long i, lx = lg(x), tx = typ(x);
    2212      146716 :   if (!is_matvec_t(tx) && tx != t_VECSMALL
    2213          49 :       && (tx != t_LIST || list_typ(x) != t_LIST_RAW)) return gcopy(x);
    2214      146695 :   if (tx == t_LIST)
    2215          28 :   { if (list_data(x)) { x = list_data(x); lx = lg(x); } else lx = 1; }
    2216      146695 :   if (lx==1) pari_err_DOMAIN("vecmin", "empty argument", "=", x,x);
    2217      146660 :   switch(typ(x))
    2218             :   {
    2219      146618 :     case t_VEC: case t_COL:
    2220      146618 :       i = vecindexmin(x); if (pi) *pi = utoipos(i);
    2221      146618 :       return gcopy(gel(x,i));
    2222          21 :     case t_MAT: {
    2223          21 :       long j, i0 = 1, j0 = 1, lx2 = lgcols(x);
    2224             :       GEN s;
    2225          21 :       if (lx2 == 1) pari_err_DOMAIN("vecmin", "empty argument", "=", x,x);
    2226          21 :       s = gcoeff(x,i0,j0); i = 2;
    2227          63 :       for (j=1; j<lx; j++,i=1)
    2228             :       {
    2229          42 :         GEN c = gel(x,j);
    2230         105 :         for (; i<lx2; i++)
    2231          63 :           if (gcmp(gel(c,i),s) < 0) { s = gel(c,i); j0=j; i0=i; }
    2232             :       }
    2233          21 :       if (pi) *pi = mkvec2(utoipos(i0), utoipos(j0));
    2234          21 :       return gcopy(s);
    2235             :     }
    2236          21 :     case t_VECSMALL:
    2237          21 :       i = vecsmall_indexmin(x); if (pi) *pi = utoipos(i);
    2238          21 :       return stoi(x[i]);
    2239             :   }
    2240             :   return NULL;/*LCOV_EXCL_LINE*/
    2241             : }
    2242             : 
    2243             : GEN
    2244       66168 : vecmax(GEN x) { return vecmax0(x, NULL); }
    2245             : GEN
    2246       34548 : vecmin(GEN x) { return vecmin0(x, NULL); }
    2247             : 
    2248             : /*******************************************************************/
    2249             : /*                                                                 */
    2250             : /*                      AFFECT long --> GEN                        */
    2251             : /*         affect long s to GEN x. Useful for initialization.      */
    2252             : /*                                                                 */
    2253             : /*******************************************************************/
    2254             : 
    2255             : static void
    2256           0 : padicaff0(GEN x)
    2257             : {
    2258           0 :   if (signe(gel(x,4)))
    2259             :   {
    2260           0 :     x[1] = evalvalp(valp(x)+precp(x));
    2261           0 :     affsi(0,gel(x,4));
    2262             :   }
    2263           0 : }
    2264             : 
    2265             : void
    2266       91840 : gaffsg(long s, GEN x)
    2267             : {
    2268       91840 :   switch(typ(x))
    2269             :   {
    2270       90930 :     case t_INT: affsi(s,x); break;
    2271         910 :     case t_REAL: affsr(s,x); break;
    2272           0 :     case t_INTMOD: modsiz(s,gel(x,1),gel(x,2)); break;
    2273           0 :     case t_FRAC: affsi(s,gel(x,1)); affsi(1,gel(x,2)); break;
    2274           0 :     case t_COMPLEX: gaffsg(s,gel(x,1)); gaffsg(0,gel(x,2)); break;
    2275           0 :     case t_PADIC: {
    2276             :       long vx;
    2277             :       GEN y;
    2278           0 :       if (!s) { padicaff0(x); break; }
    2279           0 :       vx = Z_pvalrem(stoi(s), gel(x,2), &y);
    2280           0 :       setvalp(x,vx); modiiz(y,gel(x,3),gel(x,4));
    2281           0 :       break;
    2282             :     }
    2283           0 :     case t_QUAD: gaffsg(s,gel(x,2)); gaffsg(0,gel(x,3)); break;
    2284           0 :     default: pari_err_TYPE2("=",stoi(s),x);
    2285             :   }
    2286       91840 : }
    2287             : 
    2288             : /*******************************************************************/
    2289             : /*                                                                 */
    2290             : /*                     GENERIC AFFECTATION                         */
    2291             : /*         Affect the content of x to y, whenever possible         */
    2292             : /*                                                                 */
    2293             : /*******************************************************************/
    2294             : /* x PADIC, Y INT, return lift(x * Mod(1,Y)) */
    2295             : GEN
    2296        4466 : padic_to_Fp(GEN x, GEN Y) {
    2297        4466 :   pari_sp av = avma;
    2298        4466 :   GEN p = gel(x,2), z;
    2299        4466 :   long vy, vx = valp(x);
    2300        4466 :   if (!signe(Y)) pari_err_INV("padic_to_Fp",Y);
    2301        4466 :   vy = Z_pvalrem(Y,p, &z);
    2302        4466 :   if (vx < 0 || !gequal1(z)) pari_err_OP("",x, mkintmod(gen_1,Y));
    2303        4445 :   if (vx >= vy) { set_avma(av); return gen_0; }
    2304        4088 :   z = gel(x,4);
    2305        4088 :   if (!signe(z) || vy > vx + precp(x)) pari_err_OP("",x, mkintmod(gen_1,Y));
    2306        4088 :   if (vx) z = mulii(z, powiu(p,vx));
    2307        4088 :   return gerepileuptoint(av, remii(z, Y));
    2308             : }
    2309             : ulong
    2310      421705 : padic_to_Fl(GEN x, ulong Y) {
    2311      421705 :   GEN p = gel(x,2);
    2312             :   ulong u, z;
    2313      421705 :   long vy, vx = valp(x);
    2314      421705 :   vy = u_pvalrem(Y,p, &u);
    2315      421706 :   if (vx < 0 || u != 1) pari_err_OP("",x, mkintmodu(1,Y));
    2316             :   /* Y = p^vy */
    2317      421706 :   if (vx >= vy) return 0;
    2318      375996 :   z = umodiu(gel(x,4), Y);
    2319      375999 :   if (!z || vy > vx + precp(x)) pari_err_OP("",x, mkintmodu(1,Y));
    2320      375998 :   if (vx) {
    2321           0 :     ulong pp = p[2];
    2322           0 :     z = Fl_mul(z, upowuu(pp,vx), Y); /* p^vx < p^vy = Y */
    2323             :   }
    2324      375998 :   return z;
    2325             : }
    2326             : 
    2327             : static void
    2328           0 : croak(const char *s) {
    2329             :   char *t;
    2330           0 :   t = stack_sprintf("gaffect [overwriting universal object: %s]",s);
    2331           0 :   pari_err_BUG(t);
    2332           0 : }
    2333             : 
    2334             : void
    2335      664578 : gaffect(GEN x, GEN y)
    2336             : {
    2337      664578 :   long vx, i, lx, ly, tx = typ(x), ty = typ(y);
    2338             :   pari_sp av;
    2339             :   GEN p1, num, den;
    2340             : 
    2341      664578 :   if (tx == ty) switch(tx) {
    2342      215629 :     case t_INT:
    2343      572738 :       if (!is_universal_constant(y)) { affii(x,y); return; }
    2344             :       /* y = gen_0, gnil, gen_1 or gen_2 */
    2345           0 :       if (y==gen_0)  croak("gen_0");
    2346           0 :       if (y==gen_1)  croak("gen_1");
    2347           0 :       if (y==gen_m1) croak("gen_m1");
    2348           0 :       if (y==gen_m2) croak("gen_m2");
    2349           0 :       if (y==gen_2)  croak("gen_2");
    2350           0 :       croak("gnil)");
    2351      188958 :     case t_REAL: affrr(x,y); return;
    2352           0 :     case t_INTMOD:
    2353           0 :       if (!dvdii(gel(x,1),gel(y,1))) pari_err_OP("",x,y);
    2354           0 :       modiiz(gel(x,2),gel(y,1),gel(y,2)); return;
    2355           0 :     case t_FRAC:
    2356           0 :       affii(gel(x,1),gel(y,1));
    2357           0 :       affii(gel(x,2),gel(y,2)); return;
    2358       95438 :     case t_COMPLEX:
    2359       95438 :       gaffect(gel(x,1),gel(y,1));
    2360       95438 :       gaffect(gel(x,2),gel(y,2)); return;
    2361           0 :     case t_PADIC:
    2362           0 :       if (!equalii(gel(x,2),gel(y,2))) pari_err_OP("",x,y);
    2363           0 :       modiiz(gel(x,4),gel(y,3),gel(y,4));
    2364           0 :       setvalp(y,valp(x)); return;
    2365           0 :     case t_QUAD:
    2366           0 :       if (! ZX_equal(gel(x,1),gel(y,1))) pari_err_OP("",x,y);
    2367           0 :       affii(gel(x,2),gel(y,2));
    2368           0 :       affii(gel(x,3),gel(y,3)); return;
    2369       72713 :     case t_VEC: case t_COL: case t_MAT:
    2370       72713 :       lx = lg(x); if (lx != lg(y)) pari_err_DIM("gaffect");
    2371      194584 :       for (i=1; i<lx; i++) gaffect(gel(x,i),gel(y,i));
    2372       72713 :       return;
    2373             :   }
    2374             : 
    2375             :   /* Various conversions. Avoid them, use specialized routines ! */
    2376             : 
    2377       91840 :   if (!is_const_t(ty)) pari_err_TYPE2("=",x,y);
    2378       91840 :   switch(tx)
    2379             :   {
    2380           0 :     case t_INT:
    2381             :       switch(ty)
    2382             :       {
    2383           0 :         case t_REAL:
    2384           0 :           affir(x,y); break;
    2385             : 
    2386           0 :         case t_INTMOD:
    2387           0 :           modiiz(x,gel(y,1),gel(y,2)); break;
    2388             : 
    2389           0 :         case t_COMPLEX:
    2390           0 :           gaffect(x,gel(y,1)); gaffsg(0,gel(y,2)); break;
    2391             : 
    2392           0 :         case t_PADIC:
    2393           0 :           if (!signe(x)) { padicaff0(y); break; }
    2394           0 :           av = avma;
    2395           0 :           setvalp(y, Z_pvalrem(x,gel(y,2),&p1));
    2396           0 :           affii(modii(p1,gel(y,3)), gel(y,4));
    2397           0 :           set_avma(av); break;
    2398             : 
    2399           0 :         case t_QUAD: gaffect(x,gel(y,2)); gaffsg(0,gel(y,3)); break;
    2400           0 :         default: pari_err_TYPE2("=",x,y);
    2401             :       }
    2402           0 :       break;
    2403             : 
    2404       91840 :     case t_REAL:
    2405             :       switch(ty)
    2406             :       {
    2407       91840 :         case t_COMPLEX: gaffect(x,gel(y,1)); gaffsg(0,gel(y,2)); break;
    2408           0 :         default: pari_err_TYPE2("=",x,y);
    2409             :       }
    2410       91840 :       break;
    2411             : 
    2412           0 :     case t_FRAC:
    2413             :       switch(ty)
    2414             :       {
    2415           0 :         case t_REAL: rdiviiz(gel(x,1),gel(x,2), y); break;
    2416           0 :         case t_INTMOD: av = avma;
    2417           0 :           p1 = Fp_inv(gel(x,2),gel(y,1));
    2418           0 :           affii(modii(mulii(gel(x,1),p1),gel(y,1)), gel(y,2));
    2419           0 :           set_avma(av); break;
    2420           0 :         case t_COMPLEX: gaffect(x,gel(y,1)); gaffsg(0,gel(y,2)); break;
    2421           0 :         case t_PADIC:
    2422           0 :           if (!signe(gel(x,1))) { padicaff0(y); break; }
    2423           0 :           num = gel(x,1);
    2424           0 :           den = gel(x,2);
    2425           0 :           av = avma; vx = Z_pvalrem(num, gel(y,2), &num);
    2426           0 :           if (!vx) vx = -Z_pvalrem(den,gel(y,2),&den);
    2427           0 :           setvalp(y,vx);
    2428           0 :           p1 = mulii(num,Fp_inv(den,gel(y,3)));
    2429           0 :           affii(modii(p1,gel(y,3)), gel(y,4)); set_avma(av); break;
    2430           0 :         case t_QUAD: gaffect(x,gel(y,2)); gaffsg(0,gel(y,3)); break;
    2431           0 :         default: pari_err_TYPE2("=",x,y);
    2432             :       }
    2433           0 :       break;
    2434             : 
    2435           0 :     case t_COMPLEX:
    2436           0 :       if (!gequal0(gel(x,2))) pari_err_TYPE2("=",x,y);
    2437           0 :       gaffect(gel(x,1), y);
    2438           0 :       break;
    2439             : 
    2440           0 :     case t_PADIC:
    2441             :       switch(ty)
    2442             :       {
    2443           0 :         case t_INTMOD:
    2444           0 :           av = avma; affii(padic_to_Fp(x, gel(y,1)), gel(y,2));
    2445           0 :           set_avma(av); break;
    2446           0 :         default: pari_err_TYPE2("=",x,y);
    2447             :       }
    2448           0 :       break;
    2449             : 
    2450           0 :     case t_QUAD:
    2451             :       switch(ty)
    2452             :       {
    2453           0 :         case t_INT: case t_INTMOD: case t_FRAC: case t_PADIC:
    2454           0 :           pari_err_TYPE2("=",x,y);
    2455             : 
    2456           0 :         case t_REAL:
    2457           0 :           av = avma; affgr(quadtofp(x,realprec(y)), y); set_avma(av); break;
    2458           0 :         case t_COMPLEX:
    2459           0 :           ly = precision(y); if (!ly) pari_err_TYPE2("=",x,y);
    2460           0 :           av = avma; gaffect(quadtofp(x,ly), y); set_avma(av); break;
    2461           0 :         default: pari_err_TYPE2("=",x,y);
    2462             :       }
    2463           0 :     default: pari_err_TYPE2("=",x,y);
    2464             :   }
    2465             : }
    2466             : 
    2467             : /*******************************************************************/
    2468             : /*                                                                 */
    2469             : /*           CONVERSION QUAD --> REAL, COMPLEX OR P-ADIC           */
    2470             : /*                                                                 */
    2471             : /*******************************************************************/
    2472             : GEN
    2473         252 : quadtofp(GEN x, long prec)
    2474             : {
    2475         252 :   GEN b, D, z, u = gel(x,2), v = gel(x,3);
    2476             :   pari_sp av;
    2477         252 :   if (prec < LOWDEFAULTPREC) prec = LOWDEFAULTPREC;
    2478         252 :   if (isintzero(v)) return cxcompotor(u, prec);
    2479         252 :   av = avma; D = quad_disc(x); b = gel(gel(x,1),3); /* 0 or -1 */
    2480             :   /* u + v (-b + sqrt(D)) / 2 */
    2481         252 :   if (!signe(b)) b = NULL;
    2482         252 :   if (b) u = gadd(gmul2n(u,1), v);
    2483         252 :   z = sqrtr_abs(itor(D, prec));
    2484         252 :   if (!b) shiftr_inplace(z, -1);
    2485         252 :   z = gmul(v, z);
    2486         252 :   if (signe(D) < 0)
    2487             :   {
    2488          35 :     z = mkcomplex(cxcompotor(u, prec), z);
    2489          35 :     if (!b) return gerepilecopy(av, z);
    2490           0 :     z = gmul2n(z, -1);
    2491             :   }
    2492             :   else
    2493             :   { /* if (b) x ~ (u + z) / 2 and quadnorm(x) ~ (u^2 - z^2) / 4
    2494             :      * else x ~ u + z and quadnorm(x) ~ u^2 - z^2 */
    2495         217 :     long s = gsigne(u);
    2496         217 :     if (s == -gsigne(v)) /* conjugate expression avoids cancellation */
    2497             :     {
    2498          14 :       z = gdiv(quadnorm(x), gsub(u, z));
    2499          14 :       if (b) shiftr_inplace(z, 1);
    2500             :     }
    2501             :     else
    2502             :     {
    2503         203 :       if (s) z = gadd(u, z);
    2504         203 :       if (b) shiftr_inplace(z, -1);
    2505             :     }
    2506             :   }
    2507         217 :   return gerepileupto(av, z);
    2508             : }
    2509             : 
    2510             : static GEN
    2511          28 : qtop(GEN x, GEN p, long d)
    2512             : {
    2513          28 :   GEN z, D, P, b, u = gel(x,2), v = gel(x,3);
    2514             :   pari_sp av;
    2515          28 :   if (gequal0(v)) return cvtop(u, p, d);
    2516          28 :   P = gel(x,1);
    2517          28 :   b = gel(P,3);
    2518          28 :   av = avma; D = quad_disc(x);
    2519          28 :   if (absequaliu(p,2)) d += 2;
    2520          28 :   z = Qp_sqrt(cvtop(D,p,d));
    2521          28 :   if (!z) pari_err_SQRTN("Qp_sqrt",D);
    2522          14 :   z = gmul2n(gsub(z, b), -1);
    2523             : 
    2524          14 :   z = gadd(u, gmul(v, z));
    2525          14 :   if (typ(z) != t_PADIC) /* t_INTMOD for t_QUAD of t_INTMODs... */
    2526           0 :     z = cvtop(z, p, d);
    2527          14 :   return gerepileupto(av, z);
    2528             : }
    2529             : static GEN
    2530          14 : ctop(GEN x, GEN p, long d)
    2531             : {
    2532          14 :   pari_sp av = avma;
    2533          14 :   GEN z, u = gel(x,1), v = gel(x,2);
    2534          14 :   if (isrationalzero(v)) return cvtop(u, p, d);
    2535          14 :   z = Qp_sqrt(cvtop(gen_m1, p, d - gvaluation(v, p))); /* = I */
    2536          14 :   if (!z) pari_err_SQRTN("Qp_sqrt",gen_m1);
    2537             : 
    2538          14 :   z = gadd(u, gmul(v, z));
    2539          14 :   if (typ(z) != t_PADIC) /* t_INTMOD for t_COMPLEX of t_INTMODs... */
    2540           0 :     z = cvtop(z, p, d);
    2541          14 :   return gerepileupto(av, z);
    2542             : }
    2543             : 
    2544             : /* cvtop2(stoi(s), y) */
    2545             : GEN
    2546         399 : cvstop2(long s, GEN y)
    2547             : {
    2548         399 :   GEN z, p = gel(y,2);
    2549         399 :   long v, d = signe(gel(y,4))? precp(y): 0;
    2550         399 :   if (!s) return zeropadic_shallow(p, d);
    2551         399 :   v = z_pvalrem(s, p, &s);
    2552         399 :   if (d <= 0) return zeropadic_shallow(p, v);
    2553         399 :   z = cgetg(5, t_PADIC);
    2554         399 :   z[1] = evalprecp(d) | evalvalp(v);
    2555         399 :   gel(z,2) = p;
    2556         399 :   gel(z,3) = gel(y,3);
    2557         399 :   gel(z,4) = modsi(s, gel(y,3)); return z;
    2558             : }
    2559             : 
    2560             : static GEN
    2561    17610548 : itop2_coprime(GEN x, GEN y, long v, long d)
    2562             : {
    2563    17610548 :   GEN z = cgetg(5, t_PADIC);
    2564    17609973 :   z[1] = evalprecp(d) | evalvalp(v);
    2565    17609868 :   gel(z,2) = gel(y,2);
    2566    17609868 :   gel(z,3) = gel(y,3);
    2567    17609868 :   gel(z,4) = modii(x, gel(y,3)); return z;
    2568             : }
    2569             : /* cvtop(x, gel(y,2), precp(y)), shallow */
    2570             : GEN
    2571    17618049 : cvtop2(GEN x, GEN y)
    2572             : {
    2573    17618049 :   GEN p = gel(y,2);
    2574    17618049 :   long v, d = signe(gel(y,4))? precp(y): 0;
    2575    17618049 :   switch(typ(x))
    2576             :   {
    2577    14819504 :     case t_INT:
    2578    14819504 :       if (!signe(x)) return zeropadic_shallow(p, d);
    2579    14819504 :       if (d <= 0) return zeropadic_shallow(p, Z_pval(x,p));
    2580    14815115 :       v = Z_pvalrem(x, p, &x); return itop2_coprime(x, y, v, d);
    2581             : 
    2582           0 :     case t_INTMOD:
    2583           0 :       v = Z_pval(gel(x,1),p); if (v > d) v = d;
    2584           0 :       return cvtop(gel(x,2), p, v);
    2585             : 
    2586     2797502 :     case t_FRAC:
    2587             :     {
    2588             :       GEN num, den;
    2589     2797502 :       if (d <= 0) return zeropadic_shallow(p, Q_pval(x,p));
    2590     2796221 :       num = gel(x,1); v = Z_pvalrem(num, p, &num);
    2591     2796215 :       den = gel(x,2); if (!v) v = -Z_pvalrem(den, p, &den);
    2592     2796217 :       if (!is_pm1(den)) num = mulii(num, Fp_inv(den, gel(y,3)));
    2593     2796213 :       return itop2_coprime(num, y, v, d);
    2594             :     }
    2595           7 :     case t_COMPLEX: return ctop(x, p, d);
    2596          28 :     case t_QUAD:    return qtop(x, p, d);
    2597        1197 :     case t_PADIC:
    2598        1197 :       if (!signe(gel(x,4))) return zeropadic_shallow(p, d);
    2599        1197 :       if (precp(x) <= d) return x;
    2600          35 :       return itop2_coprime(gel(x,4), y, valp(x), d); /* reduce accuracy */
    2601             :   }
    2602           0 :   pari_err_TYPE("cvtop2",x);
    2603             :   return NULL; /* LCOV_EXCL_LINE */
    2604             : }
    2605             : 
    2606             : /* assume is_const_t(tx) */
    2607             : GEN
    2608      600600 : cvtop(GEN x, GEN p, long d)
    2609             : {
    2610             :   GEN z;
    2611             :   long v;
    2612             : 
    2613      600600 :   if (typ(p) != t_INT) pari_err_TYPE("cvtop",p);
    2614      600600 :   switch(typ(x))
    2615             :   {
    2616      268324 :     case t_INT:
    2617      268324 :       if (!signe(x)) return zeropadic(p, d);
    2618      267106 :       if (d <= 0) return zeropadic(p, Z_pval(x,p));
    2619      267050 :       v = Z_pvalrem(x, p, &x);
    2620      267049 :       z = cgetg(5, t_PADIC);
    2621      267048 :       z[1] = evalprecp(d) | evalvalp(v);
    2622      267048 :       gel(z,2) = icopy(p);
    2623      267048 :       gel(z,3) = powiu(p, d);
    2624      267050 :       gel(z,4) = modii(x, gel(z,3)); return z; /* not memory-clean */
    2625             : 
    2626          28 :     case t_INTMOD:
    2627          28 :       v = Z_pval(gel(x,1),p); if (v > d) v = d;
    2628          28 :       return cvtop(gel(x,2), p, v);
    2629             : 
    2630      163019 :     case t_FRAC:
    2631             :     {
    2632             :       GEN num, den;
    2633      163019 :       if (d <= 0) return zeropadic(p, Q_pval(x,p));
    2634      163005 :       num = gel(x,1); v = Z_pvalrem(num, p, &num);
    2635      163005 :       den = gel(x,2); if (!v) v = -Z_pvalrem(den, p, &den);
    2636      163005 :       z = cgetg(5, t_PADIC);
    2637      163005 :       z[1] = evalprecp(d) | evalvalp(v);
    2638      163005 :       gel(z,2) = icopy(p);
    2639      163005 :       gel(z,3) = powiu(p, d);
    2640      163005 :       if (!is_pm1(den)) num = mulii(num, Fp_inv(den, gel(z,3)));
    2641      163004 :       gel(z,4) = modii(num, gel(z,3)); return z; /* not memory-clean */
    2642             :     }
    2643           7 :     case t_COMPLEX: return ctop(x, p, d);
    2644      169222 :     case t_PADIC:
    2645      169222 :       p = gel(x,2); /* override */
    2646      169222 :       if (!signe(gel(x,4))) return zeropadic(p, d);
    2647      169138 :       z = cgetg(5,t_PADIC);
    2648      169138 :       z[1] = x[1]; setprecp(z,d);
    2649      169138 :       gel(z,2) = icopy(p);
    2650      169138 :       gel(z,3) = powiu(p, d);
    2651      169138 :       gel(z,4) = modii(gel(x,4), gel(z,3)); return z;
    2652             : 
    2653           0 :     case t_QUAD: return qtop(x, p, d);
    2654             :   }
    2655           0 :   pari_err_TYPE("cvtop",x);
    2656             :   return NULL; /* LCOV_EXCL_LINE */
    2657             : }
    2658             : 
    2659             : GEN
    2660         126 : gcvtop(GEN x, GEN p, long r)
    2661             : {
    2662         126 :   switch(typ(x))
    2663             :   {
    2664          63 :     case t_POL: pari_APPLY_pol_normalized(gcvtop(gel(x,i),p,r));
    2665          35 :     case t_SER: pari_APPLY_ser_normalized(gcvtop(gel(x,i),p,r));
    2666           0 :     case t_POLMOD: case t_RFRAC: case t_VEC: case t_COL: case t_MAT:
    2667           0 :       pari_APPLY_same(gcvtop(gel(x,i),p,r));
    2668             :   }
    2669          98 :   return cvtop(x,p,r);
    2670             : }
    2671             : 
    2672             : long
    2673   820892523 : gexpo_safe(GEN x)
    2674             : {
    2675   820892523 :   long tx = typ(x), lx, e, f, i;
    2676             : 
    2677   820892523 :   switch(tx)
    2678             :   {
    2679   167336114 :     case t_INT:
    2680   167336114 :       return expi(x);
    2681             : 
    2682     1069840 :     case t_FRAC:
    2683     1069840 :       return expi(gel(x,1)) - expi(gel(x,2));
    2684             : 
    2685   457056235 :     case t_REAL:
    2686   457056235 :       return expo(x);
    2687             : 
    2688    87130705 :     case t_COMPLEX:
    2689    87130705 :       e = gexpo(gel(x,1));
    2690    87131139 :       f = gexpo(gel(x,2)); return maxss(e, f);
    2691             : 
    2692          91 :     case t_QUAD: {
    2693          91 :       GEN p = gel(x,1); /* mod = X^2 + {0,1}* X - {D/4, (1-D)/4})*/
    2694          91 :       long d = 1 + expi(gel(p,2))/2; /* ~ expo(sqrt(D)) */
    2695          91 :       e = gexpo(gel(x,2));
    2696          91 :       f = gexpo(gel(x,3)) + d; return maxss(e, f);
    2697             :     }
    2698    76810326 :     case t_POL: case t_SER:
    2699    76810326 :       lx = lg(x); f = -(long)HIGHEXPOBIT;
    2700   307080354 :       for (i=2; i<lx; i++) { e=gexpo(gel(x,i)); if (e>f) f=e; }
    2701    76806718 :       return f;
    2702    31607660 :     case t_VEC: case t_COL: case t_MAT:
    2703    31607660 :       lx = lg(x); f = -(long)HIGHEXPOBIT;
    2704   217109423 :       for (i=1; i<lx; i++) { e=gexpo(gel(x,i)); if (e>f) f=e; }
    2705    31607656 :       return f;
    2706             :   }
    2707          48 :   return -1-(long)HIGHEXPOBIT;
    2708             : }
    2709             : long
    2710   820650028 : gexpo(GEN x)
    2711             : {
    2712   820650028 :   long e = gexpo_safe(x);
    2713   820660613 :   if (e < -(long)HIGHEXPOBIT) pari_err_TYPE("gexpo",x);
    2714   820663495 :   return e;
    2715             : }
    2716             : GEN
    2717       89795 : gpexponent(GEN x)
    2718             : {
    2719       89795 :   long e = gexpo(x);
    2720       89795 :   return e == -(long)HIGHEXPOBIT? mkmoo(): stoi(e);
    2721             : }
    2722             : 
    2723             : long
    2724           7 : sizedigit(GEN x)
    2725             : {
    2726           7 :   return gequal0(x)? 0: (long) ((gexpo(x)+1) * LOG10_2) + 1;
    2727             : }
    2728             : 
    2729             : /* normalize series. avma is not updated */
    2730             : GEN
    2731    13353219 : normalizeser(GEN x)
    2732             : {
    2733    13353219 :   long i, lx = lg(x), vx=varn(x), vp=valser(x);
    2734             :   GEN y, z;
    2735             : 
    2736    13353219 :   if (lx == 2) { setsigne(x,0); return x; }
    2737    13352855 :   if (lx == 3) {
    2738      191358 :     z = gel(x,2);
    2739      191358 :     if (!gequal0(z)) { setsigne(x,1); return x; }
    2740       23758 :     if (isrationalzero(z)) return zeroser(vx,vp+1);
    2741        4291 :     if (isexactzero(z)) {
    2742             :       /* dangerous case: already normalized ? */
    2743         252 :       if (!signe(x)) return x;
    2744          35 :       setvalser(x,vp+1); /* no: normalize */
    2745             :     }
    2746        4074 :     setsigne(x,0); return x;
    2747             :   }
    2748    13455613 :   for (i=2; i<lx; i++)
    2749    13408881 :     if (! isrationalzero(gel(x,i))) break;
    2750    13161497 :   if (i == lx) return zeroser(vx,lx-2+vp);
    2751    13114765 :   z = gel(x,i);
    2752    13118566 :   while (i<lx && isexactzero(gel(x,i))) i++;
    2753    13114765 :   if (i == lx)
    2754             :   {
    2755         273 :     i -= 3; y = x + i;
    2756         273 :     stackdummy((pari_sp)y, (pari_sp)x);
    2757         273 :     gel(y,2) = z;
    2758         273 :     y[1] = evalsigne(0) | evalvalser(lx-2+vp) | evalvarn(vx);
    2759         273 :     y[0] = evaltyp(t_SER) | _evallg(3);
    2760         273 :     return y;
    2761             :   }
    2762             : 
    2763    13114492 :   i -= 2; y = x + i; lx -= i;
    2764    13114492 :   y[1] = evalsigne(1) | evalvalser(vp+i) | evalvarn(vx);
    2765    13114492 :   y[0] = evaltyp(t_SER) | _evallg(lx);
    2766             : 
    2767    13114492 :   stackdummy((pari_sp)y, (pari_sp)x);
    2768    13143519 :   for (i = 2; i < lx; i++)
    2769    13142630 :     if (!gequal0(gel(y, i))) return y;
    2770         889 :   setsigne(y, 0); return y;
    2771             : }
    2772             : 
    2773             : GEN
    2774           0 : normalizepol_approx(GEN x, long lx)
    2775             : {
    2776             :   long i;
    2777           0 :   for (i = lx-1; i>1; i--)
    2778           0 :     if (! gequal0(gel(x,i))) break;
    2779           0 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
    2780           0 :   setlg(x, i+1); setsigne(x, i!=1); return x;
    2781             : }
    2782             : 
    2783             : GEN
    2784   784133612 : normalizepol_lg(GEN x, long lx)
    2785             : {
    2786   784133612 :   long i, LX = 0;
    2787   784133612 :   GEN KEEP = NULL;
    2788             : 
    2789  1037965168 :   for (i = lx-1; i>1; i--)
    2790             :   {
    2791   933075081 :     GEN z = gel(x,i);
    2792   933075081 :     if (! gequal0(z) ) {
    2793   679462415 :       if (!LX) LX = i+1;
    2794   679462415 :       stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + LX));
    2795   679450442 :       x[0] = evaltyp(t_POL) | _evallg(LX);
    2796   679450442 :       setsigne(x,1); return x;
    2797   253589727 :     } else if (!isexactzero(z)) {
    2798      969328 :       if (!LX) LX = i+1; /* to be kept as leading coeff */
    2799   252867924 :     } else if (!isrationalzero(z))
    2800      832720 :       KEEP = z; /* to be kept iff all other coeffs are exact 0s */
    2801             :   }
    2802   104890087 :   if (!LX) {
    2803   104373029 :     if (KEEP) { /* e.g. Pol(Mod(0,2)) */
    2804      348358 :       gel(x,2) = KEEP;
    2805      348358 :       LX = 3;
    2806             :     } else
    2807   104024671 :       LX = 2; /* Pol(0) */
    2808             :   }
    2809   104890087 :   stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + LX));
    2810   104772496 :   x[0] = evaltyp(t_POL) | _evallg(LX);
    2811   104772496 :   setsigne(x,0); return x;
    2812             : }
    2813             : 
    2814             : /* normalize polynomial x in place */
    2815             : GEN
    2816    99198781 : normalizepol(GEN x)
    2817             : {
    2818    99198781 :   return normalizepol_lg(x, lg(x));
    2819             : }
    2820             : 
    2821             : int
    2822    79451293 : gsigne(GEN x)
    2823             : {
    2824    79451293 :   switch(typ(x))
    2825             :   {
    2826    79070970 :     case t_INT: case t_REAL: return signe(x);
    2827      379698 :     case t_FRAC: return signe(gel(x,1));
    2828         623 :     case t_QUAD:
    2829             :     {
    2830         623 :       pari_sp av = avma;
    2831         623 :       GEN T = gel(x,1), a = gel(x,2), b = gel(x,3);
    2832             :       long sa, sb;
    2833         623 :       if (signe(gel(T,2)) > 0) break;
    2834         609 :       a = gmul2n(a,1);
    2835         609 :       if (signe(gel(T,3))) a = gadd(a,b);
    2836             :       /* a + b sqrt(D) > 0 ? */
    2837         609 :       sa = gsigne(a);
    2838         609 :       sb = gsigne(b); if (sa == sb) return gc_int(av,sa);
    2839         224 :       if (sa == 0) return gc_int(av,sb);
    2840         217 :       if (sb == 0) return gc_int(av,sa);
    2841             :       /* different signs, take conjugate expression */
    2842         210 :       sb = gsigne(gsub(gsqr(a), gmul(quad_disc(x), gsqr(b))));
    2843         210 :       return gc_int(av, sb*sa);
    2844             :     }
    2845          14 :     case t_INFINITY: return inf_get_sign(x);
    2846             :   }
    2847          12 :   pari_err_TYPE("gsigne",x);
    2848             :   return 0; /* LCOV_EXCL_LINE */
    2849             : }
    2850             : 
    2851             : /*******************************************************************/
    2852             : /*                                                                 */
    2853             : /*                              LISTS                              */
    2854             : /*                                                                 */
    2855             : /*******************************************************************/
    2856             : /* make sure L can hold l elements, at least doubling the previous max number
    2857             :  * of components. */
    2858             : static void
    2859      810971 : ensure_nb(GEN L, long l)
    2860             : {
    2861      810971 :   long nmax = list_nmax(L), i, lw;
    2862             :   GEN v, w;
    2863      810971 :   if (l <= nmax) return;
    2864        1246 :   if (nmax)
    2865             :   {
    2866         490 :     nmax <<= 1;
    2867         490 :     if (l > nmax) nmax = l;
    2868         490 :     w = list_data(L); lw = lg(w);
    2869         490 :     v = newblock(nmax+1);
    2870         490 :     v[0] = w[0];
    2871     1070958 :     for (i=1; i < lw; i++) gel(v,i) = gel(w, i);
    2872         490 :     killblock(w);
    2873             :   }
    2874             :   else /* unallocated */
    2875             :   {
    2876         756 :     nmax = 32;
    2877         756 :     if (list_data(L))
    2878           0 :       pari_err(e_MISC, "store list in variable before appending elements");
    2879         756 :     v = newblock(nmax+1);
    2880         756 :     v[0] = evaltyp(t_VEC) | _evallg(1);
    2881             :   }
    2882        1246 :   list_data(L) = v;
    2883        1246 :   L[1] = evaltyp(list_typ(L))|evallg(nmax);
    2884             : }
    2885             : 
    2886             : GEN
    2887        6891 : mklist_typ(long t)
    2888             : {
    2889        6891 :   GEN L = cgetg(3,t_LIST);
    2890        6891 :   L[1] = evaltyp(t);
    2891        6891 :   list_data(L) = NULL; return L;
    2892             : }
    2893             : 
    2894             : GEN
    2895        6835 : mklist(void)
    2896             : {
    2897        6835 :   return mklist_typ(t_LIST_RAW);
    2898             : }
    2899             : 
    2900             : GEN
    2901          49 : mkmap(void)
    2902             : {
    2903          49 :   return mklist_typ(t_LIST_MAP);
    2904             : }
    2905             : 
    2906             : /* return a list with single element x, allocated on stack */
    2907             : GEN
    2908          63 : mklistcopy(GEN x)
    2909             : {
    2910          63 :   GEN y = mklist();
    2911          63 :   list_data(y) = mkveccopy(x);
    2912          63 :   return y;
    2913             : }
    2914             : 
    2915             : GEN
    2916      776125 : listput(GEN L, GEN x, long index)
    2917             : {
    2918             :   long l;
    2919             :   GEN z;
    2920             : 
    2921      776125 :   if (index < 0) pari_err_COMPONENT("listput", "<", gen_0, stoi(index));
    2922      776118 :   z = list_data(L);
    2923      776118 :   l = z? lg(z): 1;
    2924             : 
    2925      776118 :   x = gclone(x);
    2926      776118 :   if (!index || index >= l)
    2927             :   {
    2928      775964 :     ensure_nb(L, l);
    2929      775964 :     z = list_data(L); /* it may change ! */
    2930      775964 :     index = l;
    2931      775964 :     l++;
    2932             :   } else
    2933         154 :     gunclone_deep( gel(z, index) );
    2934      776118 :   gel(z,index) = x;
    2935      776118 :   z[0] = evaltyp(t_VEC) | evallg(l); /*must be after gel(z,index) is set*/
    2936      776118 :   return gel(z,index);
    2937             : }
    2938             : 
    2939             : GEN
    2940      724955 : listput0(GEN L, GEN x, long index)
    2941             : {
    2942      724955 :   if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
    2943          14 :     pari_err_TYPE("listput",L);
    2944      724941 :   (void) listput(L, x, index);
    2945      724934 :   return x;
    2946             : }
    2947             : 
    2948             : GEN
    2949       35014 : listinsert(GEN L, GEN x, long index)
    2950             : {
    2951             :   long l, i;
    2952             :   GEN z;
    2953             : 
    2954       35014 :   z = list_data(L); l = z? lg(z): 1;
    2955       35014 :   if (index <= 0) pari_err_COMPONENT("listinsert", "<=", gen_0, stoi(index));
    2956       35007 :   if (index > l) index = l;
    2957       35007 :   ensure_nb(L, l);
    2958       35007 :   BLOCK_SIGINT_START
    2959       35007 :   z = list_data(L);
    2960    87552507 :   for (i=l; i > index; i--) gel(z,i) = gel(z,i-1);
    2961       35007 :   z[0] = evaltyp(t_VEC) | evallg(l+1);
    2962       35007 :   gel(z,index) = gclone(x);
    2963       35007 :   BLOCK_SIGINT_END
    2964       35007 :   return gel(z,index);
    2965             : }
    2966             : 
    2967             : GEN
    2968       35028 : listinsert0(GEN L, GEN x, long index)
    2969             : {
    2970       35028 :   if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
    2971          14 :     pari_err_TYPE("listinsert",L);
    2972       35014 :   (void) listinsert(L, x, index);
    2973       35007 :   return x;
    2974             : }
    2975             : 
    2976             : void
    2977       21917 : listpop(GEN L, long index)
    2978             : {
    2979             :   long l, i;
    2980             :   GEN z;
    2981             : 
    2982       21917 :   if (typ(L) != t_LIST) pari_err_TYPE("listinsert",L);
    2983       21917 :   if (index < 0) pari_err_COMPONENT("listpop", "<", gen_0, stoi(index));
    2984       21917 :   z = list_data(L);
    2985       21917 :   if (!z || (l = lg(z)-1) == 0) return;
    2986             : 
    2987       21903 :   if (!index || index > l) index = l;
    2988       21903 :   BLOCK_SIGINT_START
    2989       21903 :   gunclone_deep( gel(z, index) );
    2990       21903 :   z[0] = evaltyp(t_VEC) | _evallg(l);
    2991       21910 :   for (i=index; i < l; i++) z[i] = z[i+1];
    2992       21903 :   BLOCK_SIGINT_END
    2993             : }
    2994             : 
    2995             : void
    2996          56 : listpop0(GEN L, long index)
    2997             : {
    2998          56 :   if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
    2999          14 :     pari_err_TYPE("listpop",L);
    3000          42 :   listpop(L, index);
    3001          42 : }
    3002             : 
    3003             : /* return a copy fully allocated on stack. gclone from changevalue is
    3004             :  * supposed to malloc() it */
    3005             : GEN
    3006        5930 : gtolist(GEN x)
    3007             : {
    3008             :   GEN y;
    3009             : 
    3010        5930 :   if (!x) return mklist();
    3011         370 :   switch(typ(x))
    3012             :   {
    3013         300 :     case t_VEC: case t_COL:
    3014         300 :       y = mklist();
    3015         300 :       if (lg(x) == 1) return y;
    3016         279 :       list_data(y) = gcopy(x);
    3017         279 :       settyp(list_data(y), t_VEC);
    3018         279 :       return y;
    3019           7 :     case t_LIST:
    3020           7 :       y = mklist();
    3021           7 :       list_data(y) = list_data(x)? gcopy(list_data(x)): NULL;
    3022           7 :       return y;
    3023          63 :     default:
    3024          63 :       return mklistcopy(x);
    3025             :   }
    3026             : }
    3027             : 
    3028             : void
    3029          21 : listsort(GEN L, long flag)
    3030             : {
    3031             :   long i, l;
    3032          21 :   pari_sp av = avma;
    3033             :   GEN perm, v, vnew;
    3034             : 
    3035          21 :   if (typ(L) != t_LIST) pari_err_TYPE("listsort",L);
    3036          21 :   v = list_data(L); l = v? lg(v): 1;
    3037          21 :   if (l < 3) return;
    3038          21 :   if (flag)
    3039             :   {
    3040             :     long lnew;
    3041          14 :     perm = gen_indexsort_uniq(L, (void*)&cmp_universal, cmp_nodata);
    3042          14 :     lnew = lg(perm); /* may have changed since 'uniq' */
    3043          14 :     vnew = cgetg(lnew,t_VEC);
    3044          56 :     for (i=1; i<lnew; i++) {
    3045          42 :       long c = perm[i];
    3046          42 :       gel(vnew,i) = gel(v,c);
    3047          42 :       gel(v,c) = NULL;
    3048             :     }
    3049          14 :     if (l != lnew) { /* was shortened */
    3050         105 :       for (i=1; i<l; i++)
    3051          91 :         if (gel(v,i)) gunclone_deep(gel(v,i));
    3052          14 :       l = lnew;
    3053             :     }
    3054             :   }
    3055             :   else
    3056             :   {
    3057           7 :     perm = gen_indexsort(L, (void*)&cmp_universal, cmp_nodata);
    3058           7 :     vnew = cgetg(l,t_VEC);
    3059          63 :     for (i=1; i<l; i++) gel(vnew,i) = gel(v,perm[i]);
    3060             :   }
    3061         119 :   for (i=1; i<l; i++) gel(v,i) = gel(vnew,i);
    3062          21 :   v[0] = vnew[0]; set_avma(av);
    3063             : }

Generated by: LCOV version 1.16