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 - kernel/none - level1.h (source / functions) Hit Total Coverage
Test: PARI/GP v2.16.1 lcov report (development 28904-c3aa21e911) Lines: 601 762 78.9 %
Date: 2023-12-04 07:51:13 Functions: 216 283 76.3 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : #line 2 "../src/kernel/none/level1.h"
       2             : /* Copyright (C) 2000  The PARI group.
       3             : 
       4             : This file is part of the PARI/GP package.
       5             : 
       6             : PARI/GP is free software; you can redistribute it and/or modify it under the
       7             : terms of the GNU General Public License as published by the Free Software
       8             : Foundation; either version 2 of the License, or (at your option) any later
       9             : version. It is distributed in the hope that it will be useful, but WITHOUT
      10             : ANY WARRANTY WHATSOEVER.
      11             : 
      12             : Check the License for details. You should have received a copy of it, along
      13             : with the package; see the file 'COPYING'. If not, write to the Free Software
      14             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      15             : 
      16             : /* This file defines "level 1" kernel functions.
      17             :  * These functions can be inline; they are also defined externally in
      18             :  * mpinl.c, which includes this file and never needs to be changed */
      19             : 
      20             : INLINE long
      21 82745212396 : evallg(long x)
      22             : {
      23 82745212396 :   if (x & ~LGBITS) pari_err_OVERFLOW("lg()");
      24 82750157644 :   return _evallg(x);
      25             : }
      26             : INLINE long
      27    36273326 : evalvalp(long x)
      28             : {
      29    36273326 :   long v = _evalvalp(x);
      30    36273326 :   if (v & ~VALPBITS) pari_err_OVERFLOW("valp()");
      31    36273396 :   return v;
      32             : }
      33             : INLINE long
      34    21316820 : evalvalser(long x)
      35             : {
      36    21316820 :   long v = _evalvalser(x);
      37    21316820 :   if (v & ~VALSERBITS) pari_err_OVERFLOW("valser()");
      38    21316820 :   return v;
      39             : }
      40             : INLINE long
      41 10933207496 : evalexpo(long x)
      42             : {
      43 10933207496 :   long v = _evalexpo(x);
      44 10933207496 :   if (v & ~EXPOBITS) pari_err_OVERFLOW("expo()");
      45 10935646614 :   return v;
      46             : }
      47             : INLINE long
      48    19354157 : evalprecp(long x)
      49             : {
      50    19354157 :   long v = _evalprecp(x);
      51    19354157 :   if (x & ~((1UL<<(BITS_IN_LONG-VALPnumBITS))-1)) pari_err_OVERFLOW("precp()");
      52    19354178 :   return v;
      53             : }
      54             : 
      55             : INLINE int
      56   161906958 : varncmp(long x, long y)
      57             : {
      58   161906958 :   if (varpriority[x] < varpriority[y]) return  1;
      59   125694384 :   if (varpriority[x] > varpriority[y]) return -1;
      60    68394325 :   return 0;
      61             : }
      62             : INLINE long
      63           0 : varnmin(long x, long y)
      64           0 : { return (varpriority[x] <= varpriority[y])? x: y; }
      65             : INLINE long
      66         203 : varnmax(long x, long y)
      67         203 : { return (varpriority[x] >= varpriority[y])? x: y; }
      68             : 
      69             : /* Inhibit some area gerepile-wise: declare it to be a non recursive
      70             :  * type, of length l. Thus gerepile won't inspect the zone, just copy it.
      71             :  * For the following situation:
      72             :  *   z = cgetg(t,a); av = avma; garbage(); ltop = avma;
      73             :  *   for (i=1; i<HUGE; i++) gel(z,i) = blah();
      74             :  *   stackdummy(av,ltop);
      75             :  * loses (av-ltop) words but save a costly gerepile. */
      76             : INLINE void
      77  3132460726 : stackdummy(pari_sp av, pari_sp ltop) {
      78  3132460726 :   long l = ((GEN)av) - ((GEN)ltop);
      79  3132460726 :   if (l > 0) {
      80  1064812692 :     GEN z = (GEN)ltop;
      81  1064812692 :     z[0] = evaltyp(t_VECSMALL) | evallg(l);
      82             : #ifdef DEBUG
      83             :     { long i; for (i = 1; i < l; i++) z[i] = 0; }
      84             : #endif
      85             :   }
      86  3133211399 : }
      87             : INLINE void
      88    89694310 : fixlg(GEN x, long ly) {
      89    89694310 :   long lx = lg(x), l = lx - ly;
      90    89694310 :   if (l > 0)
      91             :   { /* stackdummy(x+lx, x+ly) */
      92    47215451 :     GEN z = x + ly;
      93    47215451 :     z[0] = evaltyp(t_VECSMALL) | evallg(l);
      94    47215496 :     setlg(x, ly);
      95             : #ifdef DEBUG
      96             :     { long i; for (i = 1; i < l; i++) z[i] = 0; }
      97             : #endif
      98             :   }
      99    89694411 : }
     100             : /* update lg(z) before affrr(y, z)  [ to cater for precision loss ]*/
     101             : INLINE void
     102    42491406 : affrr_fixlg(GEN y, GEN z) { fixlg(z, lg(y)); affrr(y, z); }
     103             : 
     104             : /*******************************************************************/
     105             : /*                                                                 */
     106             : /*                       ALLOCATE ON STACK                         */
     107             : /*                                                                 */
     108             : /*******************************************************************/
     109             : INLINE ulong
     110           0 : get_avma(void) { return avma; }
     111             : INLINE void
     112 >10998*10^7 : set_avma(ulong av) { avma = av; }
     113             : 
     114             : INLINE double
     115   173122022 : gc_double(pari_sp av, double d) { set_avma(av); return d; }
     116             : INLINE long
     117   219713846 : gc_long(pari_sp av, long s) { set_avma(av); return s; }
     118             : INLINE ulong
     119    29002014 : gc_ulong(pari_sp av, ulong s) { set_avma(av); return s; }
     120             : INLINE int
     121    47535912 : gc_bool(pari_sp av, int s) { set_avma(av); return s; }
     122             : INLINE int
     123     2526127 : gc_int(pari_sp av, int s) { set_avma(av); return s; }
     124             : INLINE GEN
     125     6947022 : gc_NULL(pari_sp av) { set_avma(av); return NULL; }
     126             : INLINE GEN
     127 13688606184 : gc_const(pari_sp av, GEN x) { set_avma(av); return x; }
     128             : INLINE GEN
     129      150826 : gc_stoi(pari_sp av, long x) { set_avma(av); return stoi(x); }
     130             : INLINE GEN
     131      422896 : gc_utoi(pari_sp av, ulong x) { set_avma(av); return utoi(x); }
     132             : INLINE GEN
     133     1115617 : gc_utoipos(pari_sp av, ulong x) { set_avma(av); return utoipos(x); }
     134             : 
     135             : INLINE GEN
     136 79342175549 : new_chunk(size_t x) /* x is a number of longs */
     137             : {
     138 79342175549 :   GEN z = ((GEN) avma) - x;
     139             :   CHECK_CTRLC
     140 79342175549 :   if (x > (avma-pari_mainstack->bot) / sizeof(long))
     141          14 :     new_chunk_resize(x);
     142 79342175535 :   set_avma((pari_sp)z);
     143             : #ifdef MEMSTEP
     144             :   if (DEBUGMEM>1 && pari_mainstack->memused != DISABLE_MEMUSED) {
     145             :     long d = (long)pari_mainstack->memused - (long)z;
     146             :     if (labs(d) > 4*MEMSTEP)
     147             :     {
     148             :       pari_mainstack->memused = (pari_sp)z;
     149             :       err_printf("...%4.0lf Mbytes used\n",
     150             :                 (pari_mainstack->top-pari_mainstack->memused)/1048576.);
     151             :     }
     152             :   }
     153             : #endif
     154 79310773169 :   return z;
     155             : }
     156             : 
     157             : INLINE char *
     158    10646354 : stack_malloc(size_t N)
     159             : {
     160    10646354 :   long n = nchar2nlong(N);
     161    10646341 :   return (char*)new_chunk(n);
     162             : }
     163             : 
     164             : INLINE char *
     165    52195741 : stack_malloc_align(size_t N, long k)
     166             : {
     167    52195741 :   ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
     168    52195741 :   if (d) (void)new_chunk(d/sizeof(long));
     169    52195742 :   if (e) N += k-e;
     170    52195742 :   return (char*) new_chunk(nchar2nlong(N));
     171             : }
     172             : 
     173             : INLINE char *
     174       91786 : stack_calloc(size_t N)
     175             : {
     176       91786 :   char *p = stack_malloc(N);
     177       91785 :   memset(p, 0, N); return p;
     178             : }
     179             : 
     180             : INLINE char *
     181         956 : stack_calloc_align(size_t N, long k)
     182             : {
     183         956 :   ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
     184         956 :   if (d) (void)new_chunk(d/sizeof(long));
     185         956 :   if (e) N += k-e;
     186         956 :   return stack_calloc(N);
     187             : }
     188             : 
     189             : /* cgetg(lg(x), typ(x)), set *lx. Implicit unsetisclone() */
     190             : INLINE GEN
     191  1086553598 : cgetg_copy(GEN x, long *plx) {
     192             :   GEN y;
     193  1086553598 :   *plx = lg(x); y = new_chunk((size_t)*plx);
     194  1086575073 :   y[0] = x[0] & (TYPBITS|LGBITS); return y;
     195             : }
     196             : INLINE GEN
     197      393522 : cgetg_block(long x, long y)
     198             : {
     199      393522 :   GEN z = newblock((size_t)x);
     200      393410 :   z[0] = CLONEBIT | evaltyp(y) | evallg(x);
     201      393351 :   return z;
     202             : }
     203             : INLINE GEN
     204 21120412051 : cgetg(long x, long y)
     205             : {
     206 21120412051 :   GEN z = new_chunk((size_t)x);
     207 21108241941 :   z[0] = evaltyp(y) | evallg(x);
     208 21093717969 :   return z;
     209             : }
     210             : INLINE GEN
     211 23583749873 : cgeti(long x)
     212             : {
     213 23583749873 :   GEN z = new_chunk((size_t)x);
     214 23535664705 :   z[0] = evaltyp(t_INT) | evallg(x);
     215 23517310934 :   return z;
     216             : }
     217             : INLINE GEN
     218 13906006204 : cgetipos(long x)
     219             : {
     220 13906006204 :   GEN z = cgeti(x);
     221 13882213120 :   z[1] = evalsigne(1) | evallgefint(x);
     222 13882213120 :   return z;
     223             : }
     224             : INLINE GEN
     225   244747198 : cgetineg(long x)
     226             : {
     227   244747198 :   GEN z = cgeti(x);
     228   244748941 :   z[1] = evalsigne(-1) | evallgefint(x);
     229   244748941 :   return z;
     230             : }
     231             : INLINE GEN
     232       39653 : cgetr_block(long x)
     233             : {
     234       39653 :   long l = nbits2lg(x);
     235       39651 :   GEN z = newblock((size_t)l);
     236       39663 :   z[0] = CLONEBIT | evaltyp(t_REAL) | evallg(l);
     237       39659 :   return z;
     238             : }
     239             : INLINE GEN
     240  1513213060 : cgetr(long x)
     241             : {
     242  1513213060 :   long l = nbits2lg(x);
     243  1513127118 :   GEN z = new_chunk((size_t)l);
     244  1512616950 :   z[0] = evaltyp(t_REAL) | evallg(l);
     245  1512306186 :   return z;
     246             : }
     247             : 
     248             : /*******************************************************************/
     249             : /*                                                                 */
     250             : /*                     COPY, NEGATION, ABSOLUTE VALUE              */
     251             : /*                                                                 */
     252             : /*******************************************************************/
     253             : /* cannot do memcpy because sometimes x and y overlap */
     254             : INLINE GEN
     255  3333077171 : leafcopy(GEN x)
     256             : {
     257  3333077171 :   long lx = lg(x);
     258  3333077171 :   GEN y = new_chunk(lx); /* can't use cgetg_copy, in case x,y overlap */
     259 17508155831 :   while (--lx > 0) y[lx] = x[lx];
     260  3332630704 :   y[0] = x[0] & (TYPBITS|LGBITS); return y;
     261             : }
     262             : INLINE GEN
     263  7395244509 : icopy(GEN x)
     264             : {
     265  7395244509 :   long i = lgefint(x), lx = i;
     266  7395244509 :   GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
     267 33303299798 :   while (--i > 0) y[i] = x[i];
     268  7395186792 :   y[0] = evaltyp(t_INT) | evallg(lx);
     269  7404351606 :   return y;
     270             : }
     271             : INLINE GEN
     272   113969008 : icopyspec(GEN x, long nx)
     273             : {
     274   113969008 :   long i = nx+2, lx = i;
     275   113969008 :   GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
     276  3057856447 :   x -= 2; while (--i >= 2) y[i] = x[i];
     277   113968435 :   y[1] = evalsigne(1) | evallgefint(lx);
     278   113968435 :   y[0] = evaltyp(t_INT) | evallg(lx);
     279   113968211 :   return y;
     280             : }
     281   844033631 : INLINE GEN rcopy(GEN x) { return leafcopy(x); }
     282         707 : INLINE GEN mpcopy(GEN x) { return leafcopy(x); }
     283             : 
     284             : INLINE GEN
     285   718896399 : mpabs(GEN x) { GEN y = leafcopy(x); setabssign(y); return y; }
     286             : INLINE GEN
     287    13410998 : mpabs_shallow(GEN x) { return signe(x) < 0? mpabs(x): x; }
     288   659895752 : INLINE GEN absi(GEN x) { return mpabs(x); }
     289    40935187 : INLINE GEN absi_shallow(GEN x) { return signe(x) < 0? negi(x): x; }
     290         140 : INLINE GEN absr(GEN x) { return mpabs(x); }
     291             : 
     292             : INLINE GEN
     293   820498277 : mpneg(GEN x) { GEN y = leafcopy(x); togglesign(y); return y; }
     294   573318985 : INLINE GEN negi(GEN x) { return mpneg(x); }
     295     2512443 : INLINE GEN negr(GEN x) { return mpneg(x); }
     296             : 
     297             : /* negate in place */
     298             : INLINE void
     299  1751651203 : togglesign(GEN x) { if (x[1] & SIGNBITS) { x[1] ^= HIGHBIT; } }
     300             : INLINE void
     301   785995267 : setabssign(GEN x) { x[1] &= ~HIGHBIT; }
     302             : /* negate in place, except universal constants */
     303             : INLINE void
     304   122190654 : togglesign_safe(GEN *px)
     305             : {
     306   122190654 :   switch(*px - gen_1) /* gen_1, gen_2, gen_m1, gen_m2 */
     307             :   {
     308     2472518 :     case 0: *px = gen_m1; break;
     309           4 :     case 3: *px = gen_m2;  break;
     310      557211 :     case 6: *px = gen_1; break;
     311           0 :     case 9: *px = gen_2;  break;
     312   119160921 :     default: togglesign(*px);
     313             :   }
     314   122192882 : }
     315             : /* setsigne(y, signe(x)) */
     316             : INLINE void
     317           0 : affectsign(GEN x, GEN y)
     318             : {
     319           0 :   y[1] = (x[1] & SIGNBITS) | (y[1] & ~SIGNBITS);
     320           0 : }
     321             : /* copies sign in place, except for universal constants */
     322             : INLINE void
     323     9994817 : affectsign_safe(GEN x, GEN *py)
     324             : {
     325     9994817 :   if (((*py)[1] ^ x[1]) & HIGHBIT) togglesign_safe(py);
     326     9994817 : }
     327             : /*******************************************************************/
     328             : /*                                                                 */
     329             : /*                     GEN -> LONG, LONG -> GEN                    */
     330             : /*                                                                 */
     331             : /*******************************************************************/
     332             : /* assume x != 0, return -x as a t_INT */
     333             : INLINE GEN
     334   243898760 : utoineg(ulong x) { GEN y = cgetineg(3); y[2] = x; return y; }
     335             : /* assume x != 0, return utoi(x) */
     336             : INLINE GEN
     337 12082194951 : utoipos(ulong x) { GEN y = cgetipos(3); y[2] = x; return y; }
     338             : INLINE GEN
     339 10215628128 : utoi(ulong x) { return x? utoipos(x): gen_0; }
     340             : INLINE GEN
     341   713358566 : stoi(long x)
     342             : {
     343   713358566 :   if (!x) return gen_0;
     344   494650908 :   return x > 0? utoipos((ulong)x): utoineg((ulong)-x);
     345             : }
     346             : 
     347             : /* x 2^BIL + y */
     348             : INLINE GEN
     349  7541996754 : uutoi(ulong x, ulong y)
     350             : {
     351             :   GEN z;
     352  7541996754 :   if (!x) return utoi(y);
     353   724598802 :   z = cgetipos(4);
     354   727089748 :   *int_W_lg(z, 1, 4) = x;
     355   727089748 :   *int_W_lg(z, 0, 4) = y; return z;
     356             : }
     357             : /* - (x 2^BIL + y) */
     358             : INLINE GEN
     359      246156 : uutoineg(ulong x, ulong y)
     360             : {
     361             :   GEN z;
     362      246156 :   if (!x) return y? utoineg(y): gen_0;
     363       10425 :   z = cgetineg(4);
     364       10425 :   *int_W_lg(z, 1, 4) = x;
     365       10425 :   *int_W_lg(z, 0, 4) = y; return z;
     366             : }
     367             : 
     368             : INLINE long
     369   427419022 : itos(GEN x)
     370             : {
     371   427419022 :   long s = signe(x);
     372             :   long u;
     373             : 
     374   427419022 :   if (!s) return 0;
     375   400440388 :   u = x[2];
     376   400440388 :   if (lgefint(x) > 3 || u < 0)
     377          27 :     pari_err_OVERFLOW("t_INT-->long assignment");
     378   400438317 :   return (s>0) ? u : -u;
     379             : }
     380             : /* as itos, but return 0 if too large. Cf is_bigint */
     381             : INLINE long
     382    17526756 : itos_or_0(GEN x) {
     383             :   long n;
     384    17526756 :   if (lgefint(x) != 3 || (n = x[2]) & HIGHBIT) return 0;
     385    15782753 :   return signe(x) > 0? n: -n;
     386             : }
     387             : INLINE ulong
     388   151992306 : itou(GEN x)
     389             : {
     390   151992306 :   switch(lgefint(x)) {
     391    10044229 :     case 2: return 0;
     392   141948178 :     case 3: return x[2];
     393           0 :     default:
     394           0 :       pari_err_OVERFLOW("t_INT-->ulong assignment");
     395             :       return 0; /* LCOV_EXCL_LINE */
     396             :   }
     397             : }
     398             : 
     399             : /* as itou, but return 0 if too large. Cf is_bigint */
     400             : INLINE ulong
     401     2469097 : itou_or_0(GEN x) {
     402     2469097 :   if (lgefint(x) != 3) return 0;
     403     2448026 :   return (ulong)x[2];
     404             : }
     405             : 
     406             : INLINE ulong
     407     5520653 : umuluu_or_0(ulong x, ulong y)
     408             : {
     409             :   ulong z;
     410             :   LOCAL_HIREMAINDER;
     411     5520653 :   z = mulll(x, y);
     412     5520653 :   return hiremainder? 0: z;
     413             : }
     414             : /* return x*y if <= n, else 0. Beware overflow */
     415             : INLINE ulong
     416     5613979 : umuluu_le(ulong x, ulong y, ulong n)
     417             : {
     418             :   ulong z;
     419             :   LOCAL_HIREMAINDER;
     420     5613979 :   z = mulll(x, y);
     421     5613979 :   return (hiremainder || z > n)? 0: z;
     422             : }
     423             : 
     424             : INLINE GEN
     425   184218694 : real_0_bit(long bitprec) { GEN x=cgetg(2, t_REAL); x[1]=evalexpo(bitprec); return x; }
     426             : INLINE GEN
     427      629488 : real_0(long prec) { return real_0_bit(-prec2nbits(prec)); }
     428             : INLINE GEN
     429     2675822 : real_1_bit(long bit) { return real_1(nbits2prec(bit)); }
     430             : INLINE GEN
     431   104795035 : real_1(long prec) {
     432   104795035 :   GEN x = cgetr(prec);
     433   104786531 :   long i, l = lg(x);
     434   104786531 :   x[1] = evalsigne(1) | _evalexpo(0);
     435   474667719 :   x[2] = (long)HIGHBIT; for (i=3; i<l; i++) x[i] = 0;
     436   104786531 :   return x;
     437             : }
     438             : INLINE GEN
     439         329 : real_m1(long prec) {
     440         329 :   GEN x = cgetr(prec);
     441         329 :   long i, l = lg(x);
     442         329 :   x[1] = evalsigne(-1) | _evalexpo(0);
     443        1617 :   x[2] = (long)HIGHBIT; for (i=3; i<l; i++) x[i] = 0;
     444         329 :   return x;
     445             : }
     446             : 
     447             : /* 2.^n */
     448             : INLINE GEN
     449      591162 : real2n(long n, long prec) { GEN z = real_1(prec); setexpo(z, n); return z; }
     450             : INLINE GEN
     451           0 : real_m2n(long n, long prec) { GEN z = real_m1(prec); setexpo(z, n); return z; }
     452             : INLINE GEN
     453   413958156 : stor(long s, long prec) { GEN z = cgetr(prec); affsr(s,z); return z; }
     454             : INLINE GEN
     455    12213273 : utor(ulong s, long prec){ GEN z = cgetr(prec); affur(s,z); return z; }
     456             : INLINE GEN
     457   650697319 : itor(GEN x, long prec) { GEN z = cgetr(prec); affir(x,z); return z; }
     458             : INLINE GEN
     459   230208381 : rtor(GEN x, long prec) { GEN z = cgetr(prec); affrr(x,z); return z; }
     460             : 
     461    22121536 : INLINE ulong int_bit(GEN x, long n)
     462             : {
     463    22121536 :   long r, q = dvmdsBIL(n, &r);
     464    22109525 :   return q < lgefint(x)-2?((ulong)*int_W(x,q) >> r) & 1UL:0;
     465             : }
     466             : 
     467             : /*******************************************************************/
     468             : /*                                                                 */
     469             : /*                           COMPARISON                            */
     470             : /*                                                                 */
     471             : /*******************************************************************/
     472             : INLINE int
     473     1289805 : cmpss(long a, long b)
     474     1289805 : { return a>b? 1: (a<b? -1: 0); }
     475             : 
     476             : INLINE int
     477  1418182584 : cmpuu(ulong a, ulong b)
     478  1418182584 : { return a>b? 1: (a<b? -1: 0); }
     479             : 
     480             : INLINE int
     481     6190797 : cmpir(GEN x, GEN y)
     482             : {
     483             :   pari_sp av;
     484             :   GEN z;
     485             : 
     486     6190797 :   if (!signe(x)) return -signe(y);
     487      477245 :   if (!signe(y))
     488             :   {
     489        2248 :     if (expo(y) >= expi(x)) return 0;
     490        2220 :     return signe(x);
     491             :   }
     492      474997 :   av=avma; z = itor(x, realprec(y)); set_avma(av);
     493      474987 :   return cmprr(z,y); /* cmprr does no memory adjustment */
     494             : }
     495             : INLINE int
     496      408913 : cmpri(GEN x, GEN y) { return -cmpir(y,x); }
     497             : INLINE int
     498      586514 : cmpsr(long x, GEN y)
     499             : {
     500             :   pari_sp av;
     501             :   GEN z;
     502             : 
     503      586514 :   if (!x) return -signe(y);
     504      586514 :   av=avma; z = stor(x, LOWDEFAULTPREC); set_avma(av);
     505      586514 :   return cmprr(z,y);
     506             : }
     507             : INLINE int
     508       40996 : cmprs(GEN x, long y) { return -cmpsr(y,x); }
     509             : /* compare x and y */
     510             : INLINE int
     511     8720969 : cmpui(ulong x, GEN y)
     512             : {
     513             :   ulong p;
     514     8720969 :   if (!x) return -signe(y);
     515     8720906 :   if (signe(y) <= 0) return 1;
     516     8609137 :   if (lgefint(y) > 3) return -1;
     517     8380939 :   p = y[2]; if (p == x) return 0;
     518     8310942 :   return p < x ? 1 : -1;
     519             : }
     520             : INLINE int
     521     8716732 : cmpiu(GEN x, ulong y) { return -cmpui(y,x); }
     522             : /* compare x and |y| */
     523             : INLINE int
     524    32013873 : abscmpui(ulong x, GEN y)
     525             : {
     526    32013873 :   long l = lgefint(y);
     527             :   ulong p;
     528             : 
     529    32013873 :   if (!x) return (l > 2)? -1: 0;
     530    32013859 :   if (l == 2) return 1;
     531    31666730 :   if (l > 3) return -1;
     532    31638667 :   p = y[2]; if (p == x) return 0;
     533    31057482 :   return p < x ? 1 : -1;
     534             : }
     535             : INLINE int
     536    32012398 : abscmpiu(GEN x, ulong y) { return -abscmpui(y,x); }
     537             : INLINE int
     538     3377446 : cmpsi(long x, GEN y)
     539             : {
     540             :   ulong p;
     541             : 
     542     3377446 :   if (!x) return -signe(y);
     543             : 
     544     3374274 :   if (x > 0)
     545             :   {
     546     3373294 :     if (signe(y)<=0) return 1;
     547     3373000 :     if (lgefint(y)>3) return -1;
     548     3357692 :     p = y[2]; if (p == (ulong)x) return 0;
     549     3287962 :     return p < (ulong)x ? 1 : -1;
     550             :   }
     551             : 
     552         980 :   if (signe(y)>=0) return -1;
     553         119 :   if (lgefint(y)>3) return 1;
     554         119 :   p = y[2]; if (p == (ulong)-x) return 0;
     555          14 :   return p < (ulong)(-x) ? -1 : 1;
     556             : }
     557             : INLINE int
     558     3352335 : cmpis(GEN x, long y) { return -cmpsi(y,x); }
     559             : INLINE int
     560     2126083 : mpcmp(GEN x, GEN y)
     561             : {
     562     2126083 :   if (typ(x)==t_INT)
     563       60304 :     return (typ(y)==t_INT) ? cmpii(x,y) : cmpir(x,y);
     564     2065779 :   return (typ(y)==t_INT) ? -cmpir(y,x) : cmprr(x,y);
     565             : }
     566             : 
     567             : /* x == y ? */
     568             : INLINE int
     569     2917094 : equalui(ulong x, GEN y)
     570             : {
     571     2917094 :   if (!x) return !signe(y);
     572     2916401 :   if (signe(y) <= 0 || lgefint(y) != 3) return 0;
     573     2904837 :   return ((ulong)y[2] == (ulong)x);
     574             : }
     575             : /* x == y ? */
     576             : INLINE int
     577      524098 : equalsi(long x, GEN y)
     578             : {
     579      524098 :   if (!x) return !signe(y);
     580      524098 :   if (x > 0)
     581             :   {
     582      521899 :     if (signe(y) <= 0 || lgefint(y) != 3) return 0;
     583      469449 :     return ((ulong)y[2] == (ulong)x);
     584             :   }
     585        2199 :   if (signe(y) >= 0 || lgefint(y) != 3) return 0;
     586        2078 :   return ((ulong)y[2] == (ulong)-x);
     587             : }
     588             : /* x == |y| ? */
     589             : INLINE int
     590    38006274 : absequalui(ulong x, GEN y)
     591             : {
     592    38006274 :   if (!x) return !signe(y);
     593    38006274 :   return (lgefint(y) == 3 && (ulong)y[2] == x);
     594             : }
     595             : INLINE int
     596    36311666 : absequaliu(GEN x, ulong y) { return absequalui(y,x); }
     597             : INLINE int
     598      523916 : equalis(GEN x, long y) { return equalsi(y,x); }
     599             : INLINE int
     600     2917095 : equaliu(GEN x, ulong y) { return equalui(y,x); }
     601             : 
     602             : /* assume x != 0, is |x| == 2^n ? */
     603             : INLINE int
     604     1234799 : absrnz_equal2n(GEN x) {
     605     1234799 :   if ((ulong)x[2]==HIGHBIT)
     606             :   {
     607       29641 :     long i, lx = lg(x);
     608       92880 :     for (i = 3; i < lx; i++)
     609       67301 :       if (x[i]) return 0;
     610       25579 :     return 1;
     611             :   }
     612     1205158 :   return 0;
     613             : }
     614             : /* assume x != 0, is |x| == 1 ? */
     615             : INLINE int
     616     3684791 : absrnz_equal1(GEN x) { return !expo(x) && absrnz_equal2n(x); }
     617             : 
     618             : INLINE long
     619  8950233178 : maxss(long x, long y) { return x>y?x:y; }
     620             : INLINE long
     621  1658360620 : minss(long x, long y) { return x<y?x:y; }
     622             : INLINE long
     623     7566352 : minuu(ulong x, ulong y) { return x<y?x:y; }
     624             : INLINE long
     625    14280174 : maxuu(ulong x, ulong y) { return x>y?x:y; }
     626             : INLINE double
     627     2906588 : maxdd(double x, double y) { return x>y?x:y; }
     628             : INLINE double
     629      249264 : mindd(double x, double y) { return x<y?x:y; }
     630             : 
     631             : /*******************************************************************/
     632             : /*                                                                 */
     633             : /*                             ADD / SUB                           */
     634             : /*                                                                 */
     635             : /*******************************************************************/
     636             : INLINE GEN
     637       25067 : subuu(ulong x, ulong y)
     638             : {
     639             :   ulong z;
     640             :   LOCAL_OVERFLOW;
     641       25067 :   z = subll(x, y);
     642       25067 :   return overflow? utoineg(-z): utoi(z);
     643             : }
     644             : INLINE GEN
     645  3064210834 : adduu(ulong x, ulong y) { ulong t = x+y; return uutoi((t < x), t); }
     646             : 
     647             : INLINE GEN
     648       25067 : addss(long x, long y)
     649             : {
     650       25067 :   if (!x) return stoi(y);
     651       25067 :   if (!y) return stoi(x);
     652       25067 :   if (x > 0) return y > 0? adduu(x,y): subuu(x, -y);
     653             : 
     654       25067 :   if (y > 0) return subuu(y, -x);
     655             :   else { /* - adduu(-x, -y) */
     656           0 :     ulong t = (-x)+(-y); return uutoineg((t < (ulong)(-x)), t);
     657             :   }
     658             : }
     659       25067 : INLINE GEN subss(long x, long y) { return addss(-y,x); }
     660             : 
     661             : INLINE GEN
     662  7413571976 : subii(GEN x, GEN y)
     663             : {
     664  7413571976 :   if (x==y) return gen_0; /* frequent with x = y = gen_0 */
     665  6089994981 :   return addii_sign(x, signe(x), y, -signe(y));
     666             : }
     667             : INLINE GEN
     668 10111530320 : addii(GEN x, GEN y) { return addii_sign(x, signe(x), y, signe(y)); }
     669             : INLINE GEN
     670  2426634914 : addrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, signe(y)); }
     671             : INLINE GEN
     672   873955177 : subrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, -signe(y)); }
     673             : INLINE GEN
     674   461619936 : addir(GEN x, GEN y) { return addir_sign(x, signe(x), y, signe(y)); }
     675             : INLINE GEN
     676     2801174 : subir(GEN x, GEN y) { return addir_sign(x, signe(x), y, -signe(y)); }
     677             : INLINE GEN
     678     5957699 : subri(GEN x, GEN y) { return addir_sign(y, -signe(y), x, signe(x)); }
     679             : INLINE GEN
     680   270932395 : addsi(long x, GEN y) { return addsi_sign(x, y, signe(y)); }
     681             : INLINE GEN
     682    91964824 : addui(ulong x, GEN y) { return addui_sign(x, y, signe(y)); }
     683             : INLINE GEN
     684     5819091 : subsi(long x, GEN y) { return addsi_sign(x, y, -signe(y)); }
     685             : INLINE GEN
     686   126536399 : subui(ulong x, GEN y) { return addui_sign(x, y, -signe(y)); }
     687             : 
     688             : /*******************************************************************/
     689             : /*                                                                 */
     690             : /*                           MOD, REM, DIV                         */
     691             : /*                                                                 */
     692             : /*******************************************************************/
     693    92282551 : INLINE ulong mod2BIL(GEN x) { return *int_LSW(x); }
     694           0 : INLINE long mod64(GEN x) { return mod2BIL(x) & 63; }
     695         259 : INLINE long mod32(GEN x) { return mod2BIL(x) & 31; }
     696      236257 : INLINE long mod16(GEN x) { return mod2BIL(x) & 15; }
     697    12779052 : INLINE long mod8(GEN x)  { return mod2BIL(x) & 7; }
     698     4083931 : INLINE long mod4(GEN x)  { return mod2BIL(x) & 3; }
     699    52828636 : INLINE long mod2(GEN x)  { return mod2BIL(x) & 1; }
     700             : INLINE int
     701    83887787 : mpodd(GEN x) { return signe(x) && mod2(x); }
     702             : /* x mod 2^n, n < BITS_IN_LONG */
     703             : INLINE ulong
     704    38538490 : umodi2n(GEN x, long n)
     705             : {
     706    38538490 :   long s = signe(x);
     707    38538490 :   const ulong _2n = 1UL << n;
     708             :   ulong m;
     709    38538490 :   if (!s) return 0;
     710    37077189 :   m = *int_LSW(x) & (_2n - 1);
     711    37077189 :   if (s < 0 && m) m = _2n - m;
     712    37077189 :   return m;
     713             : }
     714           0 : INLINE ulong Mod64(GEN x){ return umodi2n(x,6); }
     715      167650 : INLINE ulong Mod32(GEN x){ return umodi2n(x,5); }
     716      244826 : INLINE ulong Mod16(GEN x){ return umodi2n(x,4); }
     717     1780948 : INLINE ulong Mod8(GEN x) { return umodi2n(x,3); }
     718    34655164 : INLINE ulong Mod4(GEN x) { return umodi2n(x,2); }
     719     1689254 : INLINE ulong Mod2(GEN x) { return umodi2n(x,1); }
     720             : 
     721             : INLINE GEN
     722    45831509 : truedivii(GEN a,GEN b) { return truedvmdii(a,b,NULL); }
     723             : INLINE GEN
     724      248133 : truedivis(GEN a, long b) { return truedvmdis(a,b,NULL); }
     725             : INLINE GEN
     726     6197775 : truedivsi(long a, GEN b) { return truedvmdsi(a,b,NULL); }
     727             : 
     728             : INLINE GEN
     729    12113747 : divii(GEN a, GEN b) { return dvmdii(a,b,NULL); }
     730             : INLINE GEN
     731  2454451550 : remii(GEN a, GEN b) { return dvmdii(a,b,ONLY_REM); }
     732             : 
     733             : INLINE GEN
     734           0 : divss(long x, long y) { return stoi(x / y); }
     735             : INLINE GEN
     736           0 : modss(long x, long y) { return utoi(smodss(x, y)); }
     737             : INLINE GEN
     738           0 : remss(long x, long y) { return stoi(x % y); }
     739             : INLINE long
     740        6414 : smodss(long x, long y)
     741             : {
     742        6414 :   long r = x%y;
     743        6414 :   return (r >= 0)? r: labs(y) + r;
     744             : }
     745             : INLINE ulong
     746   716025521 : umodsu(long x, ulong y)
     747             : {
     748   716025521 :   return x>=0 ? x%y: Fl_neg((-x)%y, y);
     749             : }
     750             : 
     751             : INLINE long
     752           0 : sdivss_rem(long x, long y, long *r)
     753             : {
     754             :   long q;
     755             :   LOCAL_HIREMAINDER;
     756           0 :   if (!y) pari_err_INV("sdivss_rem",gen_0);
     757           0 :   hiremainder = 0; q = divll((ulong)labs(x),(ulong)labs(y));
     758           0 :   if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
     759           0 :   if (y < 0) q = -q;
     760           0 :   *r = hiremainder; return q;
     761             : }
     762             : INLINE GEN
     763           0 : divss_rem(long x, long y, long *r) { return stoi(sdivss_rem(x,y,r)); }
     764             : INLINE ulong
     765   158497970 : udivuu_rem(ulong x, ulong y, ulong *r)
     766             : {
     767   158497970 :   if (!y) pari_err_INV("udivuu_rem",gen_0);
     768   158497970 :   *r = x % y; return x / y;
     769             : }
     770             : INLINE ulong
     771     2012318 : ceildivuu(ulong a, ulong b)
     772             : {
     773             :   ulong c;
     774     2012318 :   if (!a) return 0;
     775     1834039 :   c = a / b; return (a % b)? c+1: c;
     776             : }
     777             : 
     778             : INLINE ulong
     779       13282 : uabsdivui_rem(ulong x, GEN y, ulong *r)
     780             : {
     781       13282 :   long q, s = signe(y);
     782             :   LOCAL_HIREMAINDER;
     783             : 
     784       13282 :   if (!s) pari_err_INV("uabsdivui_rem",gen_0);
     785       13282 :   if (!x || lgefint(y)>3) { *r = x; return 0; }
     786       12618 :   hiremainder=0; q = (long)divll(x, (ulong)y[2]);
     787       12618 :   if (s < 0) q = -q;
     788       12618 :   *r = hiremainder; return q;
     789             : }
     790             : 
     791             : /* assume d != 0 and |n| / d can be represented as an ulong.
     792             :  * Return |n|/d, set *r = |n| % d */
     793             : INLINE ulong
     794     8349446 : uabsdiviu_rem(GEN n, ulong d, ulong *r)
     795             : {
     796     8349446 :   switch(lgefint(n))
     797             :   {
     798           0 :     case 2: *r = 0; return 0;
     799     8349446 :     case 3:
     800             :     {
     801     8349446 :       ulong nn = n[2];
     802     8349446 :       *r = nn % d; return nn / d;
     803             :     }
     804           0 :     default: /* 4 */
     805             :     {
     806             :       ulong n1, n0, q;
     807             :       LOCAL_HIREMAINDER;
     808           0 :       n0 = *int_W(n,0);
     809           0 :       n1 = *int_W(n,1);
     810           0 :       hiremainder = n1;
     811           0 :       q = divll(n0, d);
     812           0 :       *r = hiremainder; return q;
     813             :     }
     814             :   }
     815             : }
     816             : 
     817             : INLINE long
     818    51288087 : sdivsi_rem(long x, GEN y, long *r)
     819             : {
     820    51288087 :   long q, s = signe(y);
     821             :   LOCAL_HIREMAINDER;
     822             : 
     823    51288087 :   if (!s) pari_err_INV("sdivsi_rem",gen_0);
     824    51288087 :   if (!x || lgefint(y)>3 || ((long)y[2]) < 0) { *r = x; return 0; }
     825    49333834 :   hiremainder=0; q = (long)divll(labs(x), (ulong)y[2]);
     826    49333834 :   if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
     827    49333834 :   if (s < 0) q = -q;
     828    49333834 :   *r = hiremainder; return q;
     829             : }
     830             : INLINE GEN
     831           0 : divsi_rem(long s, GEN y, long *r) { return stoi(sdivsi_rem(s,y,r)); }
     832             : 
     833             : INLINE long
     834      100943 : sdivsi(long x, GEN y)
     835             : {
     836      100943 :   long q, s = signe(y);
     837             : 
     838      100943 :   if (!s) pari_err_INV("sdivsi",gen_0);
     839      100943 :   if (!x || lgefint(y)>3 || ((long)y[2]) < 0) return 0;
     840      100838 :   q = labs(x) / y[2];
     841      100838 :   if (x < 0) q = -q;
     842      100838 :   if (s < 0) q = -q;
     843      100838 :   return q;
     844             : }
     845             : 
     846             : INLINE GEN
     847           0 : dvmdss(long x, long y, GEN *z)
     848             : {
     849             :   long r;
     850           0 :   GEN q = divss_rem(x,y, &r);
     851           0 :   *z = stoi(r); return q;
     852             : }
     853             : INLINE long
     854  5961345732 : dvmdsBIL(long n, long *r) { *r = remsBIL(n); return divsBIL(n); }
     855             : INLINE ulong
     856   164569083 : dvmduBIL(ulong n, ulong *r) { *r = remsBIL(n); return divsBIL(n); }
     857             : INLINE GEN
     858           0 : dvmdsi(long x, GEN y, GEN *z)
     859             : {
     860             :   long r;
     861           0 :   GEN q = divsi_rem(x,y, &r);
     862           0 :   *z = stoi(r); return q;
     863             : }
     864             : INLINE GEN
     865           0 : dvmdis(GEN x, long y, GEN *z)
     866             : {
     867             :   long r;
     868           0 :   GEN q = divis_rem(x,y, &r);
     869           0 :   *z = stoi(r); return q;
     870             : }
     871             : 
     872             : INLINE long
     873    21006743 : smodis(GEN x, long y)
     874             : {
     875    21006743 :   pari_sp av = avma;
     876    21006743 :   long r; (void)divis_rem(x,y, &r);
     877    21006743 :   return gc_long(av, (r >= 0)? r: labs(y) + r);
     878             : }
     879             : INLINE GEN
     880    19602368 : modis(GEN x, long y) { return stoi(smodis(x,y)); }
     881             : INLINE GEN
     882    45090067 : modsi(long x, GEN y) {
     883    45090067 :   long r; (void)sdivsi_rem(x, y, &r);
     884    45090067 :   return (r >= 0)? stoi(r): addsi_sign(r, y, 1);
     885             : }
     886             : 
     887             : INLINE ulong
     888     1511727 : umodui(ulong x, GEN y)
     889             : {
     890     1511727 :   if (!signe(y)) pari_err_INV("umodui",gen_0);
     891     1511727 :   if (!x || lgefint(y) > 3) return x;
     892      415241 :   return x % (ulong)y[2];
     893             : }
     894             : 
     895             : INLINE ulong
     896      210587 : ugcdiu(GEN x, ulong y) { return ugcd(umodiu(x,y), y); }
     897             : INLINE ulong
     898        2737 : ugcdui(ulong y, GEN x) { return ugcd(umodiu(x,y), y); }
     899             : 
     900             : INLINE GEN
     901           0 : remsi(long x, GEN y)
     902           0 : { long r; (void)sdivsi_rem(x,y, &r); return stoi(r); }
     903             : INLINE GEN
     904           0 : remis(GEN x, long y)
     905             : {
     906           0 :   pari_sp av = avma;
     907             :   long r;
     908           0 :   (void)divis_rem(x,y, &r); return gc_stoi(av, r);
     909             : }
     910             : 
     911             : INLINE GEN
     912           0 : rdivis(GEN x, long y, long prec)
     913             : {
     914           0 :   GEN z = cgetr(prec);
     915           0 :   pari_sp av = avma;
     916           0 :   affrr(divrs(itor(x,prec), y),z);
     917           0 :   set_avma(av); return z;
     918             : }
     919             : INLINE GEN
     920           0 : rdivsi(long x, GEN y, long prec)
     921             : {
     922           0 :   GEN z = cgetr(prec);
     923           0 :   pari_sp av = avma;
     924           0 :   affrr(divsr(x, itor(y,prec)), z);
     925           0 :   set_avma(av); return z;
     926             : }
     927             : INLINE GEN
     928      839647 : rdivss(long x, long y, long prec)
     929             : {
     930      839647 :   GEN z = cgetr(prec);
     931      839647 :   pari_sp av = avma;
     932      839647 :   affrr(divrs(stor(x, prec), y), z);
     933      839647 :   set_avma(av); return z;
     934             : }
     935             : 
     936             : INLINE void
     937     7801768 : rdiviiz(GEN x, GEN y, GEN z)
     938             : {
     939     7801768 :   long lz = lg(z), lx = lgefint(x), ly = lgefint(y);
     940     7801768 :   if (lx == 2) { affur(0, z); return; }
     941     7801768 :   if (ly == 3)
     942             :   {
     943     2212764 :     affir(x, z); if (signe(y) < 0) togglesign(z);
     944     2212734 :     affrr(divru(z, y[2]), z);
     945             :   }
     946     5589004 :   else if (lx > lz + 1 || ly > lz + 1)
     947             :   {
     948     1851873 :     affir(x,z); affrr(divri(z, y), z);
     949             :   }
     950             :   else
     951             :   {
     952     3737131 :     long b = prec2nbits(lg2prec(lz)) + expi(y) - expi(x) + 1;
     953     3737191 :     GEN q = divii(b > 0? shifti(x, b): x, y);
     954     3737198 :     affir(q, z); if (b > 0) shiftr_inplace(z, -b);
     955             :   }
     956     7801889 :   set_avma((ulong)z);
     957             : }
     958             : INLINE GEN
     959     7762239 : rdivii(GEN x, GEN y, long prec)
     960     7762239 : { GEN z = cgetr(prec); rdiviiz(x, y, z); return z; }
     961             : INLINE GEN
     962     7346555 : fractor(GEN x, long prec)
     963     7346555 : { return rdivii(gel(x,1), gel(x,2), prec); }
     964             : 
     965             : INLINE int
     966    16089264 : dvdii(GEN x, GEN y)
     967             : {
     968    16089264 :   pari_sp av = avma;
     969             :   GEN r;
     970    16089264 :   if (!signe(x)) return 1;
     971    14795784 :   if (!signe(y)) return 0;
     972    14795784 :   r = remii(x,y);
     973    14804888 :   return gc_bool(av, r == gen_0);
     974             : }
     975             : INLINE int
     976         371 : dvdsi(long x, GEN y)
     977             : {
     978         371 :   if (x == 0) return 1;
     979         266 :   if (!signe(y)) return 0;
     980         266 :   if (lgefint(y) != 3) return 0;
     981         259 :   return x % y[2] == 0;
     982             : }
     983             : INLINE int
     984      167195 : dvdui(ulong x, GEN y)
     985             : {
     986      167195 :   if (x == 0) return 1;
     987      167195 :   if (!signe(y)) return 0;
     988      167195 :   if (lgefint(y) != 3) return 0;
     989      156574 :   return x % y[2] == 0;
     990             : }
     991             : INLINE int
     992       33492 : dvdis(GEN x, long y)
     993       33492 : { return y? smodis(x, y) == 0: signe(x) == 0; }
     994             : INLINE int
     995      576505 : dvdiu(GEN x, ulong y)
     996      576505 : { return y? umodiu(x, y) == 0: signe(x) == 0; }
     997             : 
     998             : INLINE int
     999           0 : dvdisz(GEN x, long y, GEN z)
    1000             : {
    1001           0 :   const pari_sp av = avma;
    1002             :   long r;
    1003           0 :   GEN p1 = divis_rem(x,y, &r);
    1004           0 :   set_avma(av); if (r) return 0;
    1005           0 :   affii(p1,z); return 1;
    1006             : }
    1007             : INLINE int
    1008           0 : dvdiuz(GEN x, ulong y, GEN z)
    1009             : {
    1010           0 :   const pari_sp av = avma;
    1011             :   ulong r;
    1012           0 :   GEN p1 = absdiviu_rem(x,y, &r);
    1013           0 :   set_avma(av); if (r) return 0;
    1014           0 :   affii(p1,z); return 1;
    1015             : }
    1016             : INLINE int
    1017        5431 : dvdiiz(GEN x, GEN y, GEN z)
    1018             : {
    1019        5431 :   const pari_sp av=avma;
    1020        5431 :   GEN p2, p1 = dvmdii(x,y,&p2);
    1021        5431 :   if (signe(p2)) return gc_bool(av,0);
    1022        4230 :   affii(p1,z); return gc_bool(av,1);
    1023             : }
    1024             : 
    1025             : INLINE ulong
    1026    81821485 : remlll_pre(ulong u2, ulong u1, ulong u0, ulong n, ulong ninv)
    1027             : {
    1028    81821485 :   u1 = remll_pre(u2, u1, n, ninv);
    1029    82807139 :   return remll_pre(u1, u0, n, ninv);
    1030             : }
    1031             : 
    1032             : INLINE ulong
    1033  1985012256 : Fl_sqr_pre(ulong a, ulong p, ulong pi)
    1034             : {
    1035             :   ulong x;
    1036             :   LOCAL_HIREMAINDER;
    1037  1985012256 :   x = mulll(a,a);
    1038  1985012256 :   return remll_pre(hiremainder, x, p, pi);
    1039             : }
    1040             : 
    1041             : INLINE ulong
    1042  3318584794 : Fl_mul_pre(ulong a, ulong b, ulong p, ulong pi)
    1043             : {
    1044             :   ulong x;
    1045             :   LOCAL_HIREMAINDER;
    1046  3318584794 :   x = mulll(a,b);
    1047  3318584794 :   return remll_pre(hiremainder, x, p, pi);
    1048             : }
    1049             : 
    1050             : INLINE ulong
    1051  7051050353 : Fl_addmul_pre(ulong y0, ulong x0, ulong x1, ulong p, ulong pi)
    1052             : {
    1053             :   ulong l0, h0;
    1054             :   LOCAL_HIREMAINDER;
    1055  7051050353 :   hiremainder = y0;
    1056  7051050353 :   l0 = addmul(x0, x1); h0 = hiremainder;
    1057  7051050353 :   return remll_pre(h0, l0, p, pi);
    1058             : }
    1059             : 
    1060             : INLINE ulong
    1061    55556211 : Fl_addmulmul_pre(ulong x0, ulong y0, ulong x1, ulong y1, ulong p, ulong pi)
    1062             : {
    1063             :   ulong l0, l1, h0, h1;
    1064             :   LOCAL_OVERFLOW;
    1065             :   LOCAL_HIREMAINDER;
    1066    55556211 :   l0 = mulll(x0, y0); h0 = hiremainder;
    1067    55556211 :   l1 = mulll(x1, y1); h1 = hiremainder;
    1068    55556211 :   l0 = addll(l0, l1); h0 = addllx(h0, h1);
    1069    55556211 :   return overflow ? remlll_pre(1, h0, l0, p, pi): remll_pre(h0, l0, p, pi);
    1070             : }
    1071             : 
    1072             : INLINE ulong
    1073      219824 : Fl_ellj_pre(ulong a4, ulong a6, ulong p, ulong pi)
    1074             : {
    1075             :   /* a43 = 4 a4^3 */
    1076      219824 :   ulong a43 = Fl_double(Fl_double(
    1077             :               Fl_mul_pre(a4, Fl_sqr_pre(a4, p, pi), p, pi), p), p);
    1078             :   /* a62 = 27 a6^2 */
    1079      219829 :   ulong a62 = Fl_mul_pre(Fl_sqr_pre(a6, p, pi), 27 % p, p, pi);
    1080      219827 :   ulong z1 = Fl_mul_pre(a43, 1728 % p, p, pi);
    1081      219826 :   ulong z2 = Fl_add(a43, a62, p);
    1082      219824 :   return Fl_div(z1, z2, p);
    1083             : }
    1084             : 
    1085             : /*******************************************************************/
    1086             : /*                                                                 */
    1087             : /*                        MP (INT OR REAL)                         */
    1088             : /*                                                                 */
    1089             : /*******************************************************************/
    1090             : INLINE GEN
    1091          49 : mptrunc(GEN x) { return typ(x)==t_INT? icopy(x): truncr(x); }
    1092             : INLINE GEN
    1093           0 : mpfloor(GEN x) { return typ(x)==t_INT? icopy(x): floorr(x); }
    1094             : INLINE GEN
    1095           0 : mpceil(GEN x) { return typ(x)==t_INT? icopy(x): ceilr(x); }
    1096             : INLINE GEN
    1097     1848110 : mpround(GEN x) { return typ(x) == t_INT? icopy(x): roundr(x); }
    1098             : 
    1099             : INLINE long
    1100    32488884 : mpexpo(GEN x) { return typ(x) == t_INT? expi(x): expo(x); }
    1101             : 
    1102             : INLINE GEN
    1103   353052048 : mpadd(GEN x, GEN y)
    1104             : {
    1105   353052048 :   if (typ(x)==t_INT)
    1106    20173776 :     return (typ(y)==t_INT) ? addii(x,y) : addir(x,y);
    1107   332878272 :   return (typ(y)==t_INT) ? addir(y,x) : addrr(x,y);
    1108             : }
    1109             : INLINE GEN
    1110   229821964 : mpsub(GEN x, GEN y)
    1111             : {
    1112   229821964 :   if (typ(x)==t_INT)
    1113      500340 :     return (typ(y)==t_INT) ? subii(x,y) : subir(x,y);
    1114   229321624 :   return (typ(y)==t_INT) ? subri(x,y) : subrr(x,y);
    1115             : }
    1116             : INLINE GEN
    1117   582049751 : mpmul(GEN x, GEN y)
    1118             : {
    1119   582049751 :   if (typ(x)==t_INT)
    1120    36644415 :     return (typ(y)==t_INT) ? mulii(x,y) : mulir(x,y);
    1121   545405336 :   return (typ(y)==t_INT) ? mulir(y,x) : mulrr(x,y);
    1122             : }
    1123             : INLINE GEN
    1124    74508277 : mpsqr(GEN x) { return (typ(x)==t_INT) ? sqri(x) : sqrr(x); }
    1125             : INLINE GEN
    1126      663939 : mpdiv(GEN x, GEN y)
    1127             : {
    1128      663939 :   if (typ(x)==t_INT)
    1129      260438 :     return (typ(y)==t_INT) ? divii(x,y) : divir(x,y);
    1130      403501 :   return (typ(y)==t_INT) ? divri(x,y) : divrr(x,y);
    1131             : }
    1132             : 
    1133             : /*******************************************************************/
    1134             : /*                                                                 */
    1135             : /*                          Z/nZ, n ULONG                          */
    1136             : /*                                                                 */
    1137             : /*******************************************************************/
    1138             : INLINE ulong
    1139   439338504 : Fl_double(ulong a, ulong p)
    1140             : {
    1141   439338504 :   ulong res = a << 1;
    1142   439338504 :   return (res >= p || res < a) ? res - p : res;
    1143             : }
    1144             : INLINE ulong
    1145    88557707 : Fl_triple(ulong a, ulong p)
    1146             : {
    1147    88557707 :   ulong res = a << 1;
    1148    88557707 :   if (res >= p || res < a) res -= p;
    1149    88557707 :   res += a;
    1150    88557707 :   return (res >= p || res < a)? res - p: res;
    1151             : }
    1152             : INLINE ulong
    1153    16894940 : Fl_halve(ulong a, ulong p)
    1154             : {
    1155             :   ulong ap, ap2;
    1156    16894940 :   if ((a&1UL)==0) return a>>1;
    1157     8501052 :   ap = a + p; ap2 = ap>>1;
    1158     8501052 :   return ap>=a ? ap2: (ap2|HIGHBIT);
    1159             : }
    1160             : 
    1161             : INLINE ulong
    1162  4046460238 : Fl_add(ulong a, ulong b, ulong p)
    1163             : {
    1164  4046460238 :   ulong res = a + b;
    1165  4046460238 :   return (res >= p || res < a) ? res - p : res;
    1166             : }
    1167             : INLINE ulong
    1168   695374637 : Fl_neg(ulong x, ulong p) { return x ? p - x: 0; }
    1169             : 
    1170             : INLINE ulong
    1171  6682540676 : Fl_sub(ulong a, ulong b, ulong p)
    1172             : {
    1173  6682540676 :   ulong res = a - b;
    1174  6682540676 :   return (res > a) ? res + p: res;
    1175             : }
    1176             : 
    1177             : /* centerlift(u mod p) */
    1178             : INLINE long
    1179     3900221 : Fl_center(ulong u, ulong p, ulong ps2) { return (long) (u > ps2)? u - p: u; }
    1180             : 
    1181             : INLINE ulong
    1182  2240664208 : Fl_mul(ulong a, ulong b, ulong p)
    1183             : {
    1184             :   ulong x;
    1185             :   LOCAL_HIREMAINDER;
    1186  2240664208 :   x = mulll(a,b);
    1187  2240664208 :   if (!hiremainder) return x % p;
    1188   349512300 :   (void)divll(x,p); return hiremainder;
    1189             : }
    1190             : INLINE ulong
    1191    94223330 : Fl_sqr(ulong a, ulong p)
    1192             : {
    1193             :   ulong x;
    1194             :   LOCAL_HIREMAINDER;
    1195    94223330 :   x = mulll(a,a);
    1196    94223330 :   if (!hiremainder) return x % p;
    1197    25790798 :   (void)divll(x,p); return hiremainder;
    1198             : }
    1199             : /* don't assume that p is prime: can't special case a = 0 */
    1200             : INLINE ulong
    1201    32219864 : Fl_div(ulong a, ulong b, ulong p)
    1202    32219864 : { return Fl_mul(a, Fl_inv(b, p), p); }
    1203             : 
    1204             : /*******************************************************************/
    1205             : /*                                                                 */
    1206             : /*        DEFINED FROM EXISTING ONE EXPLOITING COMMUTATIVITY       */
    1207             : /*                                                                 */
    1208             : /*******************************************************************/
    1209             : INLINE GEN
    1210     1098741 : addri(GEN x, GEN y) { return addir(y,x); }
    1211             : INLINE GEN
    1212   145810793 : addis(GEN x, long s) { return addsi(s,x); }
    1213             : INLINE GEN
    1214    88684645 : addiu(GEN x, ulong s) { return addui(s,x); }
    1215             : INLINE GEN
    1216    10865024 : addrs(GEN x, long s) { return addsr(s,x); }
    1217             : 
    1218             : INLINE GEN
    1219   122351342 : subiu(GEN x, long y) { GEN z = subui(y, x); togglesign(z); return z; }
    1220             : INLINE GEN
    1221      169856 : subis(GEN x, long y) { return addsi(-y,x); }
    1222             : INLINE GEN
    1223    13238383 : subrs(GEN x, long y) { return addsr(-y,x); }
    1224             : 
    1225             : INLINE GEN
    1226   482365713 : mulis(GEN x, long s) { return mulsi(s,x); }
    1227             : INLINE GEN
    1228   363797973 : muliu(GEN x, ulong s) { return mului(s,x); }
    1229             : INLINE GEN
    1230     2680296 : mulru(GEN x, ulong s) { return mulur(s,x); }
    1231             : INLINE GEN
    1232    34152896 : mulri(GEN x, GEN s) { return mulir(s,x); }
    1233             : INLINE GEN
    1234     7114837 : mulrs(GEN x, long s) { return mulsr(s,x); }
    1235             : 
    1236             : /*******************************************************************/
    1237             : /*                                                                 */
    1238             : /*                  VALUATION, EXPONENT, SHIFTS                    */
    1239             : /*                                                                 */
    1240             : /*******************************************************************/
    1241             : INLINE long
    1242   174747242 : vali(GEN x)
    1243             : {
    1244             :   long i;
    1245             :   GEN xp;
    1246             : 
    1247   174747242 :   if (!signe(x)) return -1;
    1248   174670449 :   xp=int_LSW(x);
    1249   179811576 :   for (i=0; !*xp; i++) xp=int_nextW(xp);
    1250   174670449 :   return vals(*xp) + i * BITS_IN_LONG;
    1251             : }
    1252             : 
    1253             : /* assume x > 0 */
    1254             : INLINE long
    1255   678813830 : expu(ulong x) { return (BITS_IN_LONG-1) - (long)bfffo(x); }
    1256             : 
    1257             : INLINE long
    1258  1909763929 : expi(GEN x)
    1259             : {
    1260  1909763929 :   const long lx=lgefint(x);
    1261  1909763929 :   return lx==2? -(long)HIGHEXPOBIT: bit_accuracy(lx)-(long)bfffo(*int_MSW(x))-1;
    1262             : }
    1263             : 
    1264             : INLINE GEN
    1265   145720723 : shiftr(GEN x, long n)
    1266             : {
    1267   145720723 :   const long e = evalexpo(expo(x)+n);
    1268   145716946 :   const GEN y = rcopy(x);
    1269             : 
    1270   145705239 :   if (e & ~EXPOBITS) pari_err_OVERFLOW("expo()");
    1271   145706977 :   y[1] = (y[1]&~EXPOBITS) | e; return y;
    1272             : }
    1273             : INLINE GEN
    1274   105739551 : mpshift(GEN x,long s) { return (typ(x)==t_INT)?shifti(x,s):shiftr(x,s); }
    1275             : 
    1276             : /* FIXME: adapt/use mpn_[lr]shift instead */
    1277             : /* z2[imin..imax] := z1[imin..imax].f shifted left sh bits
    1278             :  * (feeding f from the right). Assume sh > 0 */
    1279             : INLINE void
    1280  6376548391 : shift_left(GEN z2, GEN z1, long imin, long imax, ulong f,  ulong sh)
    1281             : {
    1282  6376548391 :   GEN sb = z1 + imin, se = z1 + imax, te = z2 + imax;
    1283  6376548391 :   ulong l, m = BITS_IN_LONG - sh, k = f >> m;
    1284 40677767605 :   while (se > sb) {
    1285 34301219214 :     l     = *se--;
    1286 34301219214 :     *te-- = (l << sh) | k;
    1287 34301219214 :     k     = l >> m;
    1288             :   }
    1289  6376548391 :   *te = (((ulong)*se) << sh) | k;
    1290  6376548391 : }
    1291             : /* z2[imin..imax] := f.z1[imin..imax-1] shifted right sh bits
    1292             :  * (feeding f from the left). Assume sh > 0 */
    1293             : INLINE void
    1294  4926002503 : shift_right(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh)
    1295             : {
    1296  4926002503 :   GEN sb = z1 + imin, se = z1 + imax, tb = z2 + imin;
    1297  4926002503 :   ulong k, l = *sb++, m = BITS_IN_LONG - sh;
    1298  4926002503 :   *tb++ = (l >> sh) | (f << m);
    1299 23733765396 :   while (sb < se) {
    1300 18807762893 :     k     = l << m;
    1301 18807762893 :     l     = *sb++;
    1302 18807762893 :     *tb++ = (l >> sh) | k;
    1303             :   }
    1304  4926002503 : }
    1305             : 
    1306             : /* Backward compatibility. Inefficient && unused */
    1307             : extern ulong hiremainder;
    1308             : INLINE ulong
    1309           0 : shiftl(ulong x, ulong y)
    1310           0 : { hiremainder = x>>(BITS_IN_LONG-y); return (x<<y); }
    1311             : 
    1312             : INLINE ulong
    1313           0 : shiftlr(ulong x, ulong y)
    1314           0 : { hiremainder = x<<(BITS_IN_LONG-y); return (x>>y); }
    1315             : 
    1316             : INLINE void
    1317   297287726 : shiftr_inplace(GEN z, long d)
    1318             : {
    1319   297287726 :   setexpo(z, expo(z)+d);
    1320   297283822 : }
    1321             : 
    1322             : /*******************************************************************/
    1323             : /*                                                                 */
    1324             : /*                           ASSIGNMENT                            */
    1325             : /*                                                                 */
    1326             : /*******************************************************************/
    1327             : INLINE void
    1328   570664299 : affii(GEN x, GEN y)
    1329             : {
    1330   570664299 :   long lx = lgefint(x);
    1331   570664299 :   if (lg(y)<lx) pari_err_OVERFLOW("t_INT-->t_INT assignment");
    1332 36521879816 :   while (--lx) y[lx] = x[lx];
    1333   570669037 : }
    1334             : INLINE void
    1335     4194890 : affsi(long s, GEN x)
    1336             : {
    1337     4194890 :   if (!s) x[1] = evalsigne(0) | evallgefint(2);
    1338             :   else
    1339             :   {
    1340     3965626 :     if (s > 0) { x[1] = evalsigne( 1) | evallgefint(3); x[2] =  s; }
    1341     1126603 :     else       { x[1] = evalsigne(-1) | evallgefint(3); x[2] = -s; }
    1342             :   }
    1343     4194890 : }
    1344             : INLINE void
    1345    44603025 : affui(ulong u, GEN x)
    1346             : {
    1347    44603025 :   if (!u) x[1] = evalsigne(0) | evallgefint(2);
    1348    44563839 :   else  { x[1] = evalsigne(1) | evallgefint(3); x[2] = u; }
    1349    44603025 : }
    1350             : 
    1351             : INLINE void
    1352   413725414 : affsr(long x, GEN y)
    1353             : {
    1354   413725414 :   long sh, i, ly = lg(y);
    1355             : 
    1356   413725414 :   if (!x)
    1357             :   {
    1358       70896 :     y[1] = evalexpo(-bit_accuracy(ly));
    1359       70896 :     return;
    1360             :   }
    1361   413654518 :   if (x < 0) {
    1362        7973 :     x = -x; sh = bfffo(x);
    1363        7973 :     y[1] = evalsigne(-1) | _evalexpo((BITS_IN_LONG-1)-sh);
    1364             :   }
    1365             :   else
    1366             :   {
    1367   413646545 :     sh = bfffo(x);
    1368   413646545 :     y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
    1369             :   }
    1370  4191771210 :   y[2] = ((ulong)x)<<sh; for (i=3; i<ly; i++) y[i]=0;
    1371             : }
    1372             : 
    1373             : INLINE void
    1374    12215314 : affur(ulong x, GEN y)
    1375             : {
    1376    12215314 :   long sh, i, ly = lg(y);
    1377             : 
    1378    12215314 :   if (!x)
    1379             :   {
    1380     1299414 :     y[1] = evalexpo(-bit_accuracy(ly));
    1381     1299414 :     return;
    1382             :   }
    1383    10915900 :   sh = bfffo(x);
    1384    10915900 :   y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
    1385    43290293 :   y[2] = x<<sh; for (i=3; i<ly; i++) y[i] = 0;
    1386             : }
    1387             : 
    1388             : INLINE void
    1389      424072 : affiz(GEN x, GEN y) { if (typ(y)==t_INT) affii(x,y); else affir(x,y); }
    1390             : INLINE void
    1391           0 : affsz(long x, GEN y) { if (typ(y)==t_INT) affsi(x,y); else affsr(x,y); }
    1392             : INLINE void
    1393      658136 : mpaff(GEN x, GEN y) { if (typ(x)==t_INT) affiz(x, y); else affrr(x,y); }
    1394             : 
    1395             : /*******************************************************************/
    1396             : /*                                                                 */
    1397             : /*                    OPERATION + ASSIGNMENT                       */
    1398             : /*                                                                 */
    1399             : /*******************************************************************/
    1400             : 
    1401           0 : INLINE void addiiz(GEN x, GEN y, GEN z)
    1402           0 : { pari_sp av = avma; affii(addii(x,y),z); set_avma(av); }
    1403           0 : INLINE void addirz(GEN x, GEN y, GEN z)
    1404           0 : { pari_sp av = avma; affrr(addir(x,y),z); set_avma(av); }
    1405           0 : INLINE void addriz(GEN x, GEN y, GEN z)
    1406           0 : { pari_sp av = avma; affrr(addri(x,y),z); set_avma(av); }
    1407     1307078 : INLINE void addrrz(GEN x, GEN y, GEN z)
    1408     1307078 : { pari_sp av = avma; affrr(addrr(x,y),z); set_avma(av); }
    1409           0 : INLINE void addsiz(long s, GEN y, GEN z)
    1410           0 : { pari_sp av = avma; affii(addsi(s,y),z); set_avma(av); }
    1411           0 : INLINE void addsrz(long s, GEN y, GEN z)
    1412           0 : { pari_sp av = avma; affrr(addsr(s,y),z); set_avma(av); }
    1413           0 : INLINE void addssz(long s, long y, GEN z)
    1414           0 : { pari_sp av = avma; affii(addss(s,y),z); set_avma(av); }
    1415             : 
    1416           0 : INLINE void diviiz(GEN x, GEN y, GEN z)
    1417           0 : { pari_sp av = avma; affii(divii(x,y),z); set_avma(av); }
    1418           0 : INLINE void divirz(GEN x, GEN y, GEN z)
    1419           0 : { pari_sp av = avma; mpaff(divir(x,y),z); set_avma(av); }
    1420           0 : INLINE void divisz(GEN x, long y, GEN z)
    1421           0 : { pari_sp av = avma; affii(divis(x,y),z); set_avma(av); }
    1422           0 : INLINE void divriz(GEN x, GEN y, GEN z)
    1423           0 : { pari_sp av = avma; affrr(divri(x,y),z); set_avma(av); }
    1424         511 : INLINE void divrrz(GEN x, GEN y, GEN z)
    1425         511 : { pari_sp av = avma; affrr(divrr(x,y),z); set_avma(av); }
    1426           0 : INLINE void divrsz(GEN y, long s, GEN z)
    1427           0 : { pari_sp av = avma; affrr(divrs(y,s),z); set_avma(av); }
    1428           0 : INLINE void divsiz(long x, GEN y, GEN z)
    1429           0 : { long junk; affsi(sdivsi_rem(x,y,&junk), z); }
    1430           0 : INLINE void divsrz(long s, GEN y, GEN z)
    1431           0 : { pari_sp av = avma; mpaff(divsr(s,y),z); set_avma(av); }
    1432           0 : INLINE void divssz(long x, long y, GEN z)
    1433           0 : { affsi(x/y, z); }
    1434             : 
    1435           0 : INLINE void modisz(GEN y, long s, GEN z)
    1436           0 : { affsi(smodis(y,s),z); }
    1437           0 : INLINE void modsiz(long s, GEN y, GEN z)
    1438           0 : { pari_sp av = avma; affii(modsi(s,y),z); set_avma(av); }
    1439           0 : INLINE void modssz(long s, long y, GEN z)
    1440           0 : { affsi(smodss(s,y),z); }
    1441             : 
    1442           0 : INLINE void mpaddz(GEN x, GEN y, GEN z)
    1443           0 : { pari_sp av = avma; mpaff(mpadd(x,y),z); set_avma(av); }
    1444           0 : INLINE void mpsubz(GEN x, GEN y, GEN z)
    1445           0 : { pari_sp av = avma; mpaff(mpsub(x,y),z); set_avma(av); }
    1446           0 : INLINE void mpmulz(GEN x, GEN y, GEN z)
    1447           0 : { pari_sp av = avma; mpaff(mpmul(x,y),z); set_avma(av); }
    1448             : 
    1449           0 : INLINE void muliiz(GEN x, GEN y, GEN z)
    1450           0 : { pari_sp av = avma; affii(mulii(x,y),z); set_avma(av); }
    1451           0 : INLINE void mulirz(GEN x, GEN y, GEN z)
    1452           0 : { pari_sp av = avma; mpaff(mulir(x,y),z); set_avma(av); }
    1453           0 : INLINE void mulriz(GEN x, GEN y, GEN z)
    1454           0 : { pari_sp av = avma; mpaff(mulri(x,y),z); set_avma(av); }
    1455      192514 : INLINE void mulrrz(GEN x, GEN y, GEN z)
    1456      192514 : { pari_sp av = avma; affrr(mulrr(x,y),z); set_avma(av); }
    1457           0 : INLINE void mulsiz(long s, GEN y, GEN z)
    1458           0 : { pari_sp av = avma; affii(mulsi(s,y),z); set_avma(av); }
    1459           0 : INLINE void mulsrz(long s, GEN y, GEN z)
    1460           0 : { pari_sp av = avma; mpaff(mulsr(s,y),z); set_avma(av); }
    1461           0 : INLINE void mulssz(long s, long y, GEN z)
    1462           0 : { pari_sp av = avma; affii(mulss(s,y),z); set_avma(av); }
    1463             : 
    1464           0 : INLINE void remiiz(GEN x, GEN y, GEN z)
    1465           0 : { pari_sp av = avma; affii(remii(x,y),z); set_avma(av); }
    1466           0 : INLINE void remisz(GEN y, long s, GEN z)
    1467           0 : { pari_sp av = avma; affii(remis(y,s),z); set_avma(av); }
    1468           0 : INLINE void remsiz(long s, GEN y, GEN z)
    1469           0 : { pari_sp av = avma; affii(remsi(s,y),z); set_avma(av); }
    1470           0 : INLINE void remssz(long s, long y, GEN z)
    1471           0 : { pari_sp av = avma; affii(remss(s,y),z); set_avma(av); }
    1472             : 
    1473          28 : INLINE void subiiz(GEN x, GEN y, GEN z)
    1474          28 : { pari_sp av = avma; affii(subii(x,y),z); set_avma(av); }
    1475           0 : INLINE void subirz(GEN x, GEN y, GEN z)
    1476           0 : { pari_sp av = avma; affrr(subir(x,y),z); set_avma(av); }
    1477           0 : INLINE void subisz(GEN y, long s, GEN z)
    1478           0 : { pari_sp av = avma; affii(addsi(-s,y),z); set_avma(av); }
    1479           0 : INLINE void subriz(GEN x, GEN y, GEN z)
    1480           0 : { pari_sp av = avma; affrr(subri(x,y),z); set_avma(av); }
    1481     1296706 : INLINE void subrrz(GEN x, GEN y, GEN z)
    1482     1296706 : { pari_sp av = avma; affrr(subrr(x,y),z); set_avma(av); }
    1483           0 : INLINE void subrsz(GEN y, long s, GEN z)
    1484           0 : { pari_sp av = avma; affrr(addsr(-s,y),z); set_avma(av); }
    1485           0 : INLINE void subsiz(long s, GEN y, GEN z)
    1486           0 : { pari_sp av = avma; affii(subsi(s,y),z); set_avma(av); }
    1487           0 : INLINE void subsrz(long s, GEN y, GEN z)
    1488           0 : { pari_sp av = avma; affrr(subsr(s,y),z); set_avma(av); }
    1489           0 : INLINE void subssz(long x, long y, GEN z) { addssz(x,-y,z); }
    1490             : 
    1491             : INLINE void
    1492           0 : dvmdssz(long x, long y, GEN z, GEN t) {
    1493           0 :   pari_sp av = avma;
    1494             :   long r;
    1495           0 :   affii(divss_rem(x,y, &r), z); set_avma(av); affsi(r,t);
    1496           0 : }
    1497             : INLINE void
    1498           0 : dvmdsiz(long x, GEN y, GEN z, GEN t) {
    1499           0 :   pari_sp av = avma;
    1500             :   long r;
    1501           0 :   affii(divsi_rem(x,y, &r), z); set_avma(av); affsi(r,t);
    1502           0 : }
    1503             : INLINE void
    1504           0 : dvmdisz(GEN x, long y, GEN z, GEN t) {
    1505           0 :   pari_sp av = avma;
    1506             :   long r;
    1507           0 :   affii(divis_rem(x,y, &r),z); set_avma(av); affsi(r,t);
    1508           0 : }
    1509             : INLINE void
    1510           0 : dvmdiiz(GEN x, GEN y, GEN z, GEN t) {
    1511           0 :   pari_sp av = avma;
    1512             :   GEN r;
    1513           0 :   affii(dvmdii(x,y,&r),z); affii(r,t); set_avma(av);
    1514           0 : }

Generated by: LCOV version 1.14