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 - language - init.c (source / functions) Hit Total Coverage
Test: PARI/GP v2.12.1 lcov report (development 25819-e703fe1174) Lines: 1100 1431 76.9 %
Date: 2020-09-18 06:10:04 Functions: 131 159 82.4 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /* Copyright (C) 2000-2003  The PARI group.
       2             : 
       3             : This file is part of the PARI/GP package.
       4             : 
       5             : PARI/GP is free software; you can redistribute it and/or modify it under the
       6             : terms of the GNU General Public License as published by the Free Software
       7             : Foundation. It is distributed in the hope that it will be useful, but WITHOUT
       8             : ANY WARRANTY WHATSOEVER.
       9             : 
      10             : Check the License for details. You should have received a copy of it, along
      11             : with the package; see the file 'COPYING'. If not, write to the Free Software
      12             : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
      13             : 
      14             : /*******************************************************************/
      15             : /*                                                                 */
      16             : /*        INITIALIZING THE SYSTEM, ERRORS, STACK MANAGEMENT        */
      17             : /*                                                                 */
      18             : /*******************************************************************/
      19             : /* _GNU_SOURCE is needed before first include to get RUSAGE_THREAD */
      20             : #undef _GNU_SOURCE /* avoid warning */
      21             : #define _GNU_SOURCE
      22             : #include <string.h>
      23             : #if defined(_WIN32) || defined(__CYGWIN32__)
      24             : #  include "../systems/mingw/mingw.h"
      25             : #  include <process.h>
      26             : #endif
      27             : #include "paricfg.h"
      28             : #if defined(STACK_CHECK) && !defined(__EMX__) && !defined(_WIN32)
      29             : #  include <sys/types.h>
      30             : #  include <sys/time.h>
      31             : #  include <sys/resource.h>
      32             : #endif
      33             : #if defined(HAS_WAITPID) && defined(HAS_SETSID)
      34             : #  include <sys/wait.h>
      35             : #endif
      36             : #ifdef HAS_MMAP
      37             : #  include <sys/mman.h>
      38             : #endif
      39             : #if defined(USE_GETTIMEOFDAY) || defined(USE_GETRUSAGE) || defined(USE_TIMES)
      40             : #  include <sys/time.h>
      41             : #endif
      42             : #if defined(USE_GETRUSAGE)
      43             : #  include <sys/resource.h>
      44             : #endif
      45             : #if defined(USE_FTIME) || defined(USE_FTIMEFORWALLTIME)
      46             : #  include <sys/timeb.h>
      47             : #endif
      48             : #if defined(USE_CLOCK_GETTIME) || defined(USE_TIMES)
      49             : #  include <time.h>
      50             : #endif
      51             : #if defined(USE_TIMES)
      52             : #  include <sys/times.h>
      53             : #endif
      54             : #include "pari.h"
      55             : #include "paripriv.h"
      56             : #include "anal.h"
      57             : 
      58             : const double LOG10_2 = 0.3010299956639812; /* log_10(2) */
      59             : const double LOG2_10 = 3.321928094887362;  /* log_2(10) */
      60             : 
      61             : GEN gnil, gen_0, gen_1, gen_m1, gen_2, gen_m2, ghalf, err_e_STACK;
      62             : 
      63             : static const ulong readonly_constants[] = {
      64             :   evaltyp(t_INT) | _evallg(2),  /* gen_0 */
      65             :   evallgefint(2),
      66             :   evaltyp(t_INT) | _evallg(2),  /* gnil */
      67             :   evallgefint(2),
      68             :   evaltyp(t_INT) | _evallg(3),  /* gen_1 */
      69             :   evalsigne(1) | evallgefint(3),
      70             :   1,
      71             :   evaltyp(t_INT) | _evallg(3),  /* gen_2 */
      72             :   evalsigne(1) | evallgefint(3),
      73             :   2,
      74             :   evaltyp(t_INT) | _evallg(3),  /* gen_m1 */
      75             :   evalsigne(-1) | evallgefint(3),
      76             :   1,
      77             :   evaltyp(t_INT) | _evallg(3),  /* gen_m2 */
      78             :   evalsigne(-1) | evallgefint(3),
      79             :   2,
      80             :   evaltyp(t_ERROR) | _evallg(2), /* err_e_STACK */
      81             :   e_STACK,
      82             :   evaltyp(t_FRAC) | _evallg(3), /* ghalf */
      83             :   (ulong)(readonly_constants+4),
      84             :   (ulong)(readonly_constants+7)
      85             : };
      86             : THREAD GEN zetazone, bernzone, primetab;
      87             : byteptr diffptr;
      88             : FILE    *pari_outfile, *pari_errfile, *pari_logfile, *pari_infile;
      89             : char    *current_logfile, *current_psfile, *pari_datadir;
      90             : long    gp_colors[c_LAST];
      91             : int     disable_color;
      92             : ulong   DEBUGFILES, DEBUGLEVEL, DEBUGMEM;
      93             : long    DEBUGVAR;
      94             : ulong   pari_mt_nbthreads;
      95             : long    precreal;
      96             : ulong   precdl, pari_logstyle;
      97             : gp_data *GP_DATA;
      98             : 
      99             : entree  **varentries;
     100             : THREAD long *varpriority;
     101             : 
     102             : THREAD pari_sp avma;
     103             : THREAD struct pari_mainstack *pari_mainstack;
     104             : 
     105             : static void ** MODULES;
     106             : static pari_stack s_MODULES;
     107             : const long functions_tblsz = 135; /* size of functions_hash */
     108             : entree **functions_hash, **defaults_hash;
     109             : 
     110             : char *(*cb_pari_fgets_interactive)(char *s, int n, FILE *f);
     111             : int (*cb_pari_get_line_interactive)(const char*, const char*, filtre_t *F);
     112             : void (*cb_pari_quit)(long);
     113             : void (*cb_pari_init_histfile)(void);
     114             : void (*cb_pari_ask_confirm)(const char *);
     115             : int  (*cb_pari_handle_exception)(long);
     116             : int  (*cb_pari_err_handle)(GEN);
     117             : int  (*cb_pari_whatnow)(PariOUT *out, const char *, int);
     118             : void (*cb_pari_sigint)(void);
     119             : void (*cb_pari_pre_recover)(long);
     120             : void (*cb_pari_err_recover)(long);
     121             : int (*cb_pari_break_loop)(int);
     122             : int (*cb_pari_is_interactive)(void);
     123             : void (*cb_pari_start_output)();
     124             : 
     125             : const char * pari_library_path = NULL;
     126             : 
     127             : static THREAD GEN global_err_data;
     128             : THREAD jmp_buf *iferr_env;
     129             : const long CATCH_ALL = -1;
     130             : 
     131             : static void pari_init_timer(void);
     132             : 
     133             : /*********************************************************************/
     134             : /*                                                                   */
     135             : /*                       BLOCKS & CLONES                             */
     136             : /*                                                                   */
     137             : /*********************************************************************/
     138             : /*#define DEBUG*/
     139             : static THREAD long next_block;
     140             : static THREAD GEN cur_block; /* current block in block list */
     141             : static THREAD GEN root_block; /* current block in block list */
     142             : #ifdef DEBUG
     143             : static THREAD long NUM;
     144             : #endif
     145             : 
     146             : static void
     147      177818 : pari_init_blocks(void)
     148             : {
     149      177818 :   next_block = 0; cur_block = NULL; root_block = NULL;
     150             : #ifdef DEBUG
     151             :   NUM = 0;
     152             : #endif
     153      177818 : }
     154             : 
     155             : static void
     156      174475 : pari_close_blocks(void)
     157             : {
     158     6766838 :   while (cur_block) killblock(cur_block);
     159      176600 : }
     160             : 
     161             : static long
     162 10139361824 : blockheight(GEN bl) { return bl? bl_height(bl): 0; }
     163             : 
     164             : static long
     165  2438929933 : blockbalance(GEN bl)
     166  2438929933 : { return bl ? blockheight(bl_left(bl)) - blockheight(bl_right(bl)): 0; }
     167             : 
     168             : static void
     169  2630760162 : fix_height(GEN bl)
     170  2630760162 : { bl_height(bl) = maxss(blockheight(bl_left(bl)), blockheight(bl_right(bl)))+1; }
     171             : 
     172             : static GEN
     173    54780824 : bl_rotright(GEN y)
     174             : {
     175    54780824 :   GEN x = bl_left(y), t = bl_right(x);
     176    54780824 :   bl_right(x) = y;
     177    54780824 :   bl_left(y)  = t;
     178    54780824 :   fix_height(y);
     179    54780583 :   fix_height(x);
     180    54780415 :   return x;
     181             : }
     182             : 
     183             : static GEN
     184    50811743 : bl_rotleft(GEN x)
     185             : {
     186    50811743 :   GEN y = bl_right(x), t = bl_left(y);
     187    50811743 :   bl_left(y)  = x;
     188    50811743 :   bl_right(x) = t;
     189    50811743 :   fix_height(x);
     190    50811682 :   fix_height(y);
     191    50811632 :   return y;
     192             : }
     193             : 
     194             : static GEN
     195  1490414012 : blockinsert(GEN x, GEN bl, long *d)
     196             : {
     197             :   long b, c;
     198  1490414012 :   if (!bl)
     199             :   {
     200   211724245 :     bl_left(x)=NULL; bl_right(x)=NULL;
     201   211724245 :     bl_height(x)=1; return x;
     202             :   }
     203  1278689767 :   c = cmpuu((ulong)x, (ulong)bl);
     204  1278691077 :   if (c < 0)
     205   533797090 :     bl_left(bl) = blockinsert(x, bl_left(bl), d);
     206   744893987 :   else if (c > 0)
     207   744893987 :     bl_right(bl) = blockinsert(x, bl_right(bl), d);
     208           0 :   else return bl; /* ??? Already exist in the tree ? */
     209  1278689832 :   fix_height(bl);
     210  1278689184 :   b = blockbalance(bl);
     211  1278689767 :   if (b > 1)
     212             :   {
     213    36020403 :     if (*d > 0) bl_left(bl) = bl_rotleft(bl_left(bl));
     214    36020403 :     return bl_rotright(bl);
     215             :   }
     216  1242669364 :   if (b < -1)
     217             :   {
     218    20858190 :     if (*d < 0) bl_right(bl) = bl_rotright(bl_right(bl));
     219    20858190 :     return bl_rotleft(bl);
     220             :   }
     221  1221811174 :   *d = c; return bl;
     222             : }
     223             : 
     224             : static GEN
     225  1352615269 : blockdelete(GEN x, GEN bl)
     226             : {
     227             :   long b;
     228  1352615269 :   if (!bl) return NULL; /* ??? Do not exist in the tree */
     229  1352615269 :   if (x < bl)
     230   482302487 :     bl_left(bl) = blockdelete(x, bl_left(bl));
     231   870312782 :   else if (x > bl)
     232   618018521 :     bl_right(bl) = blockdelete(x, bl_right(bl));
     233             :   else
     234             :   {
     235   252294261 :     if (!bl_left(bl) && !bl_right(bl)) return NULL;
     236    89618013 :     else if (!bl_left(bl)) return bl_right(bl);
     237    60622545 :     else if (!bl_right(bl)) return bl_left(bl);
     238             :     else
     239             :     {
     240    40571619 :       GEN r = bl_right(bl);
     241    57518003 :       while (bl_left(r)) r = bl_left(r);
     242    40571619 :       bl_right(r) = blockdelete(r, bl_right(bl));
     243    40576225 :       bl_left(r) = bl_left(bl);
     244    40576225 :       bl = r;
     245             :     }
     246             :   }
     247  1140890099 :   fix_height(bl);
     248  1140886298 :   b = blockbalance(bl);
     249  1140887033 :   if (b > 1)
     250             :   {
     251    14408859 :     if (blockbalance(bl_left(bl)) >= 0) return bl_rotright(bl);
     252             :     else
     253     3445429 :     { bl_left(bl) = bl_rotleft(bl_left(bl)); return bl_rotright(bl); }
     254             :   }
     255  1126478174 :   if (b < -1)
     256             :   {
     257     4946804 :     if (blockbalance(bl_right(bl)) <= 0) return bl_rotleft(bl);
     258             :     else
     259     1123107 :     { bl_right(bl) = bl_rotright(bl_right(bl)); return bl_rotleft(bl); }
     260             :   }
     261  1121531370 :   return bl;
     262             : }
     263             : 
     264             : static GEN
     265   623343460 : blocksearch(GEN x, GEN bl)
     266             : {
     267   623343460 :   if (isclone(x)) return x;
     268   451263689 :   if (isonstack(x) || is_universal_constant(x)) return NULL;
     269   728331408 :   while (bl)
     270             :   {
     271   726053649 :     if (x >= bl  && x < bl + bl_size(bl))
     272   193416230 :       return bl;
     273   532637419 :     bl = x < bl ? bl_left(bl): bl_right(bl);
     274             :   }
     275     2277759 :   return NULL; /* Unknown address */
     276             : }
     277             : 
     278             : void
     279   314613137 : clone_lock(GEN x)
     280             : {
     281   314613137 :   GEN y = blocksearch(x, root_block);
     282   314246109 :   if (y && isclone(y))
     283             :   {
     284   185123622 :     if (DEBUGMEM > 2)
     285           0 :       err_printf("locking block no %ld: %08lx from %08lx\n", bl_num(y), y, x);
     286   185123622 :     ++bl_refc(y);
     287             :   }
     288   314246109 : }
     289             : 
     290             : void
     291   250123714 : clone_unlock(GEN x)
     292             : {
     293   250123714 :   GEN y = blocksearch(x, root_block);
     294   250029308 :   if (y && isclone(y))
     295             :   {
     296   128523357 :     if (DEBUGMEM > 2)
     297           0 :       err_printf("unlocking block no %ld: %08lx from %08lx\n", bl_num(y), y, x);
     298   128523357 :     gunclone(y);
     299             :   }
     300   250029308 : }
     301             : 
     302             : void
     303    58976165 : clone_unlock_deep(GEN x)
     304             : {
     305    58976165 :   GEN y = blocksearch(x, root_block);
     306    58976165 :   if (y && isclone(y))
     307             :   {
     308    51848724 :     if (DEBUGMEM > 2)
     309           0 :       err_printf("unlocking deep block no %ld: %08lx from %08lx\n", bl_num(y), y, x);
     310    51848724 :     gunclone_deep(y);
     311             :   }
     312    58976165 : }
     313             : 
     314             : /* Return x, where:
     315             :  * x[-8]: AVL height
     316             :  * x[-7]: adress of left child or NULL
     317             :  * x[-6]: adress of right child or NULL
     318             :  * x[-5]: size
     319             :  * x[-4]: reference count
     320             :  * x[-3]: adress of next block
     321             :  * x[-2]: adress of preceding block.
     322             :  * x[-1]: number of allocated blocs.
     323             :  * x[0..n-1]: malloc-ed memory. */
     324             : GEN
     325   211723782 : newblock(size_t n)
     326             : {
     327   211723782 :   long d = 0;
     328   211723782 :   long *x = (long *) pari_malloc((n + BL_HEAD)*sizeof(long)) + BL_HEAD;
     329             : 
     330   211726145 :   bl_size(x) = n;
     331   211726145 :   bl_refc(x) = 1;
     332   211726145 :   bl_next(x) = NULL;
     333   211726145 :   bl_prev(x) = cur_block;
     334   211726145 :   bl_num(x)  = next_block++;
     335   211726145 :   if (cur_block) bl_next(cur_block) = x;
     336   211726145 :   root_block = blockinsert(x, root_block, &d);
     337             : #ifdef DEBUG
     338             :   err_printf("+ %ld\n", ++NUM);
     339             : #endif
     340   211723513 :   if (DEBUGMEM > 2)
     341           0 :     err_printf("new block, size %6lu (no %ld): %08lx\n", n, next_block-1, x);
     342   211723667 :   return cur_block = x;
     343             : }
     344             : 
     345             : GEN
     346       31259 : gcloneref(GEN x)
     347             : {
     348       31259 :   if (isclone(x)) { ++bl_refc(x); return x; }
     349       30692 :   else return gclone(x);
     350             : }
     351             : 
     352             : void
     353           0 : gclone_refc(GEN x) { ++bl_refc(x); }
     354             : 
     355             : void
     356   345436135 : gunclone(GEN x)
     357             : {
     358   345436135 :   if (--bl_refc(x) > 0) return;
     359   211717986 :   BLOCK_SIGINT_START;
     360   211723674 :   root_block = blockdelete(x, root_block);
     361   211714825 :   if (bl_next(x)) bl_prev(bl_next(x)) = bl_prev(x);
     362             :   else
     363             :   {
     364    29489057 :     cur_block = bl_prev(x);
     365    29489057 :     next_block = bl_num(x);
     366             :   }
     367   211714825 :   if (bl_prev(x)) bl_next(bl_prev(x)) = bl_next(x);
     368   211714825 :   if (DEBUGMEM > 2)
     369           0 :     err_printf("killing block (no %ld): %08lx\n", bl_num(x), x);
     370   211714825 :   free((void*)bl_base(x)); /* pari_free not needed: we already block */
     371   211714825 :   BLOCK_SIGINT_END;
     372             : #ifdef DEBUG
     373             :   err_printf("- %ld\n", NUM--);
     374             : #endif
     375             : }
     376             : 
     377             : /* Recursively look for clones in the container and kill them. Then kill
     378             :  * container if clone. SIGINT could be blocked until it returns */
     379             : void
     380  2931200376 : gunclone_deep(GEN x)
     381             : {
     382             :   long i, lx;
     383             :   GEN v;
     384  2931200376 :   if (isclone(x) && bl_refc(x) > 1) { --bl_refc(x); return; }
     385  2879349626 :   BLOCK_SIGINT_START;
     386  2879349628 :   switch(typ(x))
     387             :   {
     388   118007203 :     case t_VEC: case t_COL: case t_MAT:
     389   118007203 :       lx = lg(x);
     390  2799463020 :       for (i=1;i<lx;i++) gunclone_deep(gel(x,i));
     391   118007141 :       break;
     392         862 :     case t_LIST:
     393         862 :       v = list_data(x); lx = v? lg(v): 1;
     394      789206 :       for (i=1;i<lx;i++) gunclone_deep(gel(v,i));
     395         862 :       if (v) killblock(v);
     396         862 :       break;
     397             :   }
     398  2879349566 :   if (isclone(x)) gunclone(x);
     399  2879349218 :   BLOCK_SIGINT_END;
     400             : }
     401             : 
     402             : int
     403      192760 : pop_entree_block(entree *ep, long loc)
     404             : {
     405      192760 :   GEN x = (GEN)ep->value;
     406      192760 :   if (bl_num(x) < loc) return 0; /* older */
     407         350 :   if (DEBUGMEM>2)
     408           0 :     err_printf("popping %s (block no %ld)\n", ep->name, bl_num(x));
     409         350 :   gunclone_deep(x); return 1;
     410             : }
     411             : 
     412             : /***************************************************************************
     413             :  **                                                                       **
     414             :  **                           Export                                      **
     415             :  **                                                                       **
     416             :  ***************************************************************************/
     417             : 
     418             : static hashtable *export_hash;
     419             : static void
     420        1690 : pari_init_export(void)
     421             : {
     422        1690 :   export_hash = hash_create_str(1,0);
     423        1690 : }
     424             : static void
     425        1680 : pari_close_export(void)
     426             : {
     427        1680 :   hash_destroy(export_hash);
     428        1680 : }
     429             : 
     430             : /* Exported values are blocks, but do not have the clone bit set so that they
     431             :  * are not affected by clone_lock and ensure_nb, etc. */
     432             : 
     433             : void
     434          52 : export_add(const char *str, GEN val)
     435             : {
     436             :   hashentry *h;
     437          52 :   val = gclone(val); unsetisclone(val);
     438          52 :   h = hash_search(export_hash, (void*) str);
     439          52 :   if (h)
     440             :   {
     441          21 :     GEN v = (GEN)h->val;
     442          21 :     h->val = val;
     443          21 :     setisclone(v); gunclone(v);
     444             :   }
     445             :   else
     446          31 :     hash_insert(export_hash,(void*)str, (void*) val);
     447          52 : }
     448             : 
     449             : void
     450          24 : export_del(const char *str)
     451             : {
     452          24 :   hashentry *h = hash_remove(export_hash,(void*)str);
     453          24 :   if (h)
     454             :   {
     455          24 :     GEN v = (GEN)h->val;
     456          24 :     setisclone(v); gunclone(v);
     457          24 :     pari_free(h);
     458             :   }
     459          24 : }
     460             : 
     461             : GEN
     462         576 : export_get(const char *str)
     463             : {
     464         576 :   return hash_haskey_GEN(export_hash,(void*)str);
     465             : }
     466             : 
     467             : void
     468           6 : unexportall(void)
     469             : {
     470           6 :   pari_sp av = avma;
     471           6 :   GEN keys = hash_keys(export_hash);
     472           6 :   long i, l = lg(keys);
     473          24 :   for (i = 1; i < l; i++) mt_export_del((const char *)keys[i]);
     474           6 :   set_avma(av);
     475           6 : }
     476             : 
     477             : void
     478           6 : exportall(void)
     479             : {
     480             :   long i;
     481         816 :   for (i = 0; i < functions_tblsz; i++)
     482             :   {
     483             :     entree *ep;
     484        8262 :     for (ep = functions_hash[i]; ep; ep = ep->next)
     485        7452 :       if (EpVALENCE(ep)==EpVAR) mt_export_add(ep->name, (GEN)ep->value);
     486             :   }
     487           6 : }
     488             : 
     489             : /*********************************************************************/
     490             : /*                                                                   */
     491             : /*                       C STACK SIZE CONTROL                        */
     492             : /*                                                                   */
     493             : /*********************************************************************/
     494             : /* Avoid core dump on deep recursion. Adapted Perl code by Dominic Dunlop */
     495             : THREAD void *PARI_stack_limit = NULL;
     496             : 
     497             : #ifdef STACK_CHECK
     498             : 
     499             : #  ifdef __EMX__                                /* Emulate */
     500             : void
     501             : pari_stackcheck_init(void *pari_stack_base)
     502             : {
     503             :   if (!pari_stack_base) { PARI_stack_limit = NULL; return; }
     504             :   PARI_stack_limit = get_stack(1./16, 32*1024);
     505             : }
     506             : #  elif _WIN32
     507             : void
     508             : pari_stackcheck_init(void *pari_stack_base)
     509             : {
     510             :   ulong size = 1UL << 21;
     511             :   if (!pari_stack_base) { PARI_stack_limit = NULL; return; }
     512             :   if (size > (ulong)pari_stack_base)
     513             :     PARI_stack_limit = (void*)(((ulong)pari_stack_base) / 16);
     514             :   else
     515             :     PARI_stack_limit = (void*)((ulong)pari_stack_base - (size/16)*15);
     516             : }
     517             : #  else /* !__EMX__ && !_WIN32 */
     518             : /* Set PARI_stack_limit to (a little above) the lowest safe address that can be
     519             :  * used on the stack. Leave PARI_stack_limit at its initial value (NULL) to
     520             :  * show no check should be made [init failed]. Assume stack grows downward. */
     521             : void
     522      179482 : pari_stackcheck_init(void *pari_stack_base)
     523             : {
     524             :   struct rlimit rip;
     525             :   ulong size;
     526      179482 :   if (!pari_stack_base) { PARI_stack_limit = NULL; return; }
     527      179482 :   if (getrlimit(RLIMIT_STACK, &rip)) return;
     528      179580 :   size = rip.rlim_cur;
     529      179580 :   if (size == (ulong)RLIM_INFINITY || size > (ulong)pari_stack_base)
     530           0 :     PARI_stack_limit = (void*)(((ulong)pari_stack_base) / 16);
     531             :   else
     532      179588 :     PARI_stack_limit = (void*)((ulong)pari_stack_base - (size/16)*15);
     533             : }
     534             : #  endif /* !__EMX__ */
     535             : 
     536             : #else
     537             : void
     538             : pari_stackcheck_init(void *pari_stack_base)
     539             : {
     540             :   (void) pari_stack_base; PARI_stack_limit = NULL;
     541             : }
     542             : #endif /* STACK_CHECK */
     543             : 
     544             : /*******************************************************************/
     545             : /*                         HEAP TRAVERSAL                          */
     546             : /*******************************************************************/
     547             : struct getheap_t { long n, l; };
     548             : /* x is a block, not necessarily a clone [x[0] may not be set] */
     549             : static void
     550        6720 : f_getheap(GEN x, void *D)
     551             : {
     552        6720 :   struct getheap_t *T = (struct getheap_t*)D;
     553        6720 :   T->n++;
     554        6720 :   T->l += bl_size(x) + BL_HEAD;
     555        6720 : }
     556             : GEN
     557          84 : getheap(void)
     558             : {
     559          84 :   struct getheap_t T = { 0, 0 };
     560          84 :   traverseheap(&f_getheap, &T); return mkvec2s(T.n, T.l);
     561             : }
     562             : 
     563             : static void
     564       13524 : traverseheap_r(GEN bl, void(*f)(GEN, void *), void *data)
     565             : {
     566       13524 :   if (!bl) return;
     567        6720 :   traverseheap_r(bl_left(bl), f, data);
     568        6720 :   traverseheap_r(bl_right(bl), f, data);
     569        6720 :   f(bl, data);
     570             : }
     571             : 
     572             : void
     573          84 : traverseheap( void(*f)(GEN, void *), void *data)
     574             : {
     575          84 :   traverseheap_r(root_block,f, data);
     576          84 : }
     577             : 
     578             : /*********************************************************************/
     579             : /*                          DAEMON / FORK                            */
     580             : /*********************************************************************/
     581             : #if defined(HAS_WAITPID) && defined(HAS_SETSID)
     582             : /* Properly fork a process, detaching from main process group without creating
     583             :  * zombies on exit. Parent returns 1, son returns 0 */
     584             : int
     585          60 : pari_daemon(void)
     586             : {
     587          60 :   pid_t pid = fork();
     588         120 :   switch(pid) {
     589           0 :       case -1: return 1; /* father, fork failed */
     590          60 :       case 0:
     591          60 :         (void)setsid(); /* son becomes process group leader */
     592          60 :         if (fork()) _exit(0); /* now son exits, also when fork fails */
     593           0 :         break; /* grandson: its father is the son, which exited,
     594             :                 * hence father becomes 'init', that'll take care of it */
     595          60 :       default: /* father, fork succeeded */
     596          60 :         (void)waitpid(pid,NULL,0); /* wait for son to exit, immediate */
     597          60 :         return 1;
     598             :   }
     599             :   /* grandson. The silly '!' avoids a gcc-8 warning (unused value) */
     600           0 :   (void)!freopen("/dev/null","r",stdin);
     601           0 :   return 0;
     602             : }
     603             : #else
     604             : int
     605             : pari_daemon(void)
     606             : {
     607             :   pari_err_IMPL("pari_daemon without waitpid & setsid");
     608             :   return 0;
     609             : }
     610             : #endif
     611             : 
     612             : /*********************************************************************/
     613             : /*                                                                   */
     614             : /*                       SYSTEM INITIALIZATION                       */
     615             : /*                                                                   */
     616             : /*********************************************************************/
     617             : static int try_to_recover = 0;
     618             : THREAD VOLATILE int PARI_SIGINT_block = 0, PARI_SIGINT_pending = 0;
     619             : 
     620             : /*********************************************************************/
     621             : /*                         SIGNAL HANDLERS                           */
     622             : /*********************************************************************/
     623             : static void
     624           0 : dflt_sigint_fun(void) { pari_err(e_MISC, "user interrupt"); }
     625             : 
     626             : #if defined(_WIN32) || defined(__CYGWIN32__)
     627             : int win32ctrlc = 0, win32alrm = 0;
     628             : void
     629             : dowin32ctrlc(void)
     630             : {
     631             :   win32ctrlc = 0;
     632             :   cb_pari_sigint();
     633             : }
     634             : #endif
     635             : 
     636             : static void
     637           0 : pari_handle_SIGINT(void)
     638             : {
     639             : #ifdef _WIN32
     640             :   if (++win32ctrlc >= 5) _exit(3);
     641             : #else
     642           0 :   cb_pari_sigint();
     643             : #endif
     644           0 : }
     645             : 
     646             : typedef void (*pari_sighandler_t)(int);
     647             : 
     648             : pari_sighandler_t
     649       18530 : os_signal(int sig, pari_sighandler_t f)
     650             : {
     651             : #ifdef HAS_SIGACTION
     652             :   struct sigaction sa, oldsa;
     653             : 
     654       18530 :   sa.sa_handler = f;
     655       18530 :   sigemptyset(&sa.sa_mask);
     656       18530 :   sa.sa_flags = SA_NODEFER;
     657             : 
     658       18530 :   if (sigaction(sig, &sa, &oldsa)) return NULL;
     659       18530 :   return oldsa.sa_handler;
     660             : #else
     661             :   return signal(sig,f);
     662             : #endif
     663             : }
     664             : 
     665             : void
     666           0 : pari_sighandler(int sig)
     667             : {
     668             :   const char *msg;
     669             : #ifndef HAS_SIGACTION
     670             :   /*SYSV reset the signal handler in the handler*/
     671             :   (void)os_signal(sig,pari_sighandler);
     672             : #endif
     673           0 :   switch(sig)
     674             :   {
     675             : #ifdef SIGBREAK
     676             :     case SIGBREAK:
     677             :       if (PARI_SIGINT_block==1)
     678             :       {
     679             :         PARI_SIGINT_pending=SIGBREAK;
     680             :         mt_sigint();
     681             :       }
     682             :       else pari_handle_SIGINT();
     683             :       return;
     684             : #endif
     685             : 
     686             : #ifdef SIGINT
     687           0 :     case SIGINT:
     688           0 :       if (PARI_SIGINT_block==1)
     689             :       {
     690           0 :         PARI_SIGINT_pending=SIGINT;
     691           0 :         mt_sigint();
     692             :       }
     693           0 :       else pari_handle_SIGINT();
     694           0 :       return;
     695             : #endif
     696             : 
     697             : #ifdef SIGSEGV
     698           0 :     case SIGSEGV:
     699           0 :       msg="PARI/GP (Segmentation Fault)"; break;
     700             : #endif
     701             : #ifdef SIGBUS
     702           0 :     case SIGBUS:
     703           0 :       msg="PARI/GP (Bus Error)"; break;
     704             : #endif
     705             : #ifdef SIGFPE
     706           0 :     case SIGFPE:
     707           0 :       msg="PARI/GP (Floating Point Exception)"; break;
     708             : #endif
     709             : 
     710             : #ifdef SIGPIPE
     711           0 :     case SIGPIPE:
     712             :     {
     713           0 :       pariFILE *f = GP_DATA->pp->file;
     714           0 :       if (f && pari_outfile == f->file)
     715             :       {
     716           0 :         GP_DATA->pp->file = NULL; /* to avoid oo recursion on error */
     717           0 :         pari_outfile = stdout; pari_fclose(f);
     718           0 :         pari_err(e_MISC, "Broken Pipe, resetting file stack...");
     719             :       }
     720             :       return; /* LCOV_EXCL_LINE */
     721             :     }
     722             : #endif
     723             : 
     724           0 :     default: msg="signal handling"; break;
     725             :   }
     726           0 :   pari_err_BUG(msg);
     727             : }
     728             : 
     729             : void
     730        3370 : pari_sig_init(void (*f)(int))
     731             : {
     732             : #ifdef SIGBUS
     733        3370 :   (void)os_signal(SIGBUS,f);
     734             : #endif
     735             : #ifdef SIGFPE
     736        3370 :   (void)os_signal(SIGFPE,f);
     737             : #endif
     738             : #ifdef SIGINT
     739        3370 :   (void)os_signal(SIGINT,f);
     740             : #endif
     741             : #ifdef SIGBREAK
     742             :   (void)os_signal(SIGBREAK,f);
     743             : #endif
     744             : #ifdef SIGPIPE
     745        3370 :   (void)os_signal(SIGPIPE,f);
     746             : #endif
     747             : #ifdef SIGSEGV
     748        3370 :   (void)os_signal(SIGSEGV,f);
     749             : #endif
     750        3370 : }
     751             : 
     752             : /*********************************************************************/
     753             : /*                      STACK AND UNIVERSAL CONSTANTS                */
     754             : /*********************************************************************/
     755             : static void
     756        1690 : init_universal_constants(void)
     757             : {
     758        1690 :   gen_0  = (GEN)readonly_constants;
     759        1690 :   gnil   = (GEN)readonly_constants+2;
     760        1690 :   gen_1  = (GEN)readonly_constants+4;
     761        1690 :   gen_2  = (GEN)readonly_constants+7;
     762        1690 :   gen_m1 = (GEN)readonly_constants+10;
     763        1690 :   gen_m2 = (GEN)readonly_constants+13;
     764        1690 :   err_e_STACK = (GEN)readonly_constants+16;
     765        1690 :   ghalf  = (GEN)readonly_constants+18;
     766        1690 : }
     767             : 
     768             : static void
     769      178048 : pari_init_errcatch(void)
     770             : {
     771      178048 :   iferr_env = NULL;
     772      178048 :   global_err_data = NULL;
     773      178048 : }
     774             : 
     775             : /*********************************************************************/
     776             : /*                           INIT DEFAULTS                           */
     777             : /*********************************************************************/
     778             : void
     779        1690 : pari_init_defaults(void)
     780             : {
     781             :   long i;
     782        1690 :   initout(1);
     783             : 
     784             : #ifdef LONG_IS_64BIT
     785        1458 :   precreal = 128;
     786             : #else
     787         232 :   precreal = 96;
     788             : #endif
     789             : 
     790        1690 :   precdl = 16;
     791        1690 :   DEBUGFILES = DEBUGLEVEL = 0;
     792        1690 :   DEBUGMEM = 1;
     793        1690 :   disable_color = 1;
     794        1690 :   pari_logstyle = logstyle_none;
     795             : 
     796        1690 :   current_psfile = pari_strdup("pari.ps");
     797        1690 :   current_logfile= pari_strdup("pari.log");
     798        1690 :   pari_logfile = NULL;
     799             : 
     800        1690 :   pari_datadir = os_getenv("GP_DATA_DIR");
     801        1690 :   if (!pari_datadir)
     802             :   {
     803             : #if defined(_WIN32) || defined(__CYGWIN32__)
     804             :     if (paricfg_datadir[0]=='@' && paricfg_datadir[1]==0)
     805             :       pari_datadir = win32_datadir();
     806             :     else
     807             : #endif
     808        1690 :       pari_datadir = pari_strdup(paricfg_datadir);
     809             :   }
     810           0 :   else pari_datadir= pari_strdup(pari_datadir);
     811       13520 :   for (i=0; i<c_LAST; i++) gp_colors[i] = c_NONE;
     812        1690 : }
     813             : 
     814             : /*********************************************************************/
     815             : /*                   FUNCTION HASHTABLES, MODULES                    */
     816             : /*********************************************************************/
     817             : extern entree functions_basic[], functions_default[];
     818             : static void
     819        1690 : pari_init_functions(void)
     820             : {
     821        1690 :   pari_stack_init(&s_MODULES, sizeof(*MODULES),(void**)&MODULES);
     822        1690 :   pari_stack_pushp(&s_MODULES,functions_basic);
     823        1690 :   functions_hash = (entree**) pari_calloc(sizeof(entree*)*functions_tblsz);
     824        1690 :   pari_fill_hashtable(functions_hash, functions_basic);
     825        1690 :   defaults_hash = (entree**) pari_calloc(sizeof(entree*)*functions_tblsz);
     826        1690 :   pari_add_defaults_module(functions_default);
     827        1690 : }
     828             : 
     829             : void
     830        1680 : pari_add_module(entree *ep)
     831             : {
     832        1680 :   pari_fill_hashtable(functions_hash, ep);
     833        1680 :   pari_stack_pushp(&s_MODULES, ep);
     834        1680 : }
     835             : 
     836             : void
     837        1690 : pari_add_defaults_module(entree *ep)
     838        1690 : { pari_fill_hashtable(defaults_hash, ep); }
     839             : 
     840             : /*********************************************************************/
     841             : /*                       PARI MAIN STACK                             */
     842             : /*********************************************************************/
     843             : 
     844             : #ifdef HAS_MMAP
     845             : #define PARI_STACK_ALIGN (sysconf(_SC_PAGE_SIZE))
     846             : #ifndef MAP_ANONYMOUS
     847             : #define MAP_ANONYMOUS MAP_ANON
     848             : #endif
     849             : #ifndef MAP_NORESERVE
     850             : #define MAP_NORESERVE 0
     851             : #endif
     852             : static void *
     853      178331 : pari_mainstack_malloc(size_t size)
     854             : {
     855             :   /* Check that the system allows reserving "size" bytes. This is just
     856             :    * a check, we immediately free the memory. */
     857      178331 :   void *b = mmap(NULL, size, PROT_READ|PROT_WRITE,
     858             :                              MAP_PRIVATE|MAP_ANONYMOUS, -1, 0);
     859      178331 :   if (b == MAP_FAILED) return NULL;
     860      178331 :   munmap(b, size);
     861             : 
     862             :   /* Map again, this time with MAP_NORESERVE. On some operating systems
     863             :    * like Cygwin, this is needed because remapping with PROT_NONE and
     864             :    * MAP_NORESERVE does not work as expected. */
     865      178331 :   b = mmap(NULL, size, PROT_READ|PROT_WRITE,
     866             :                        MAP_PRIVATE|MAP_ANONYMOUS|MAP_NORESERVE, -1, 0);
     867      178331 :   if (b == MAP_FAILED) return NULL;
     868      178331 :   return b;
     869             : }
     870             : 
     871             : static void
     872      178321 : pari_mainstack_mfree(void *s, size_t size)
     873             : {
     874      178321 :   munmap(s, size);
     875      178321 : }
     876             : 
     877             : /* Completely discard the memory mapped between the addresses "from"
     878             :  * and "to" (which must be page-aligned).
     879             :  *
     880             :  * We use mmap() with PROT_NONE, which means that the underlying memory
     881             :  * is freed and that the kernel should not commit memory for it. We
     882             :  * still keep the mapping such that we can change the flags to
     883             :  * PROT_READ|PROT_WRITE later.
     884             :  *
     885             :  * NOTE: remapping with MAP_FIXED and PROT_NONE is not the same as
     886             :  * calling mprotect(..., PROT_NONE) because the latter will keep the
     887             :  * memory committed (this is in particular relevant on Linux with
     888             :  * vm.overcommit = 2). This remains true even when calling
     889             :  * madvise(..., MADV_DONTNEED). */
     890             : static void
     891      268399 : pari_mainstack_mreset(pari_sp from, pari_sp to)
     892             : {
     893      268399 :   size_t s = to - from;
     894             :   void *addr, *res;
     895      268399 :   if (!s) return;
     896             : 
     897           0 :   addr = (void*)from;
     898           0 :   res = mmap(addr, s, PROT_NONE,
     899             :              MAP_FIXED|MAP_PRIVATE|MAP_ANONYMOUS|MAP_NORESERVE, -1, 0);
     900           0 :   if (res != addr) pari_err(e_MEM);
     901             : }
     902             : 
     903             : /* Commit (make available) the virtual memory mapped between the
     904             :  * addresses "from" and "to" (which must be page-aligned).
     905             :  * Return 0 if successful, -1 if failed. */
     906             : static int
     907      268399 : pari_mainstack_mextend(pari_sp from, pari_sp to)
     908             : {
     909      268399 :   size_t s = to - from;
     910      268399 :   return mprotect((void*)from, s, PROT_READ|PROT_WRITE);
     911             : }
     912             : 
     913             : /* Set actual stack size to the given size. This sets st->size and
     914             :  * st->bot. If not enough system memory is available, this can fail.
     915             :  * Return 1 if successful, 0 if failed (in that case, st->size is not
     916             :  * changed) */
     917             : static int
     918      268399 : pari_mainstack_setsize(struct pari_mainstack *st, size_t size)
     919             : {
     920      268399 :   pari_sp newbot = st->top - size;
     921             :   /* Align newbot to pagesize */
     922      268399 :   pari_sp alignbot = newbot & ~(pari_sp)(PARI_STACK_ALIGN - 1);
     923      268399 :   if (pari_mainstack_mextend(alignbot, st->top))
     924             :   {
     925             :     /* Making the memory available did not work: limit vsize to the
     926             :      * current actual stack size. */
     927           0 :     st->vsize = st->size;
     928           0 :     pari_warn(warnstack, st->vsize);
     929           0 :     return 0;
     930             :   }
     931      268399 :   pari_mainstack_mreset(st->vbot, alignbot);
     932      268399 :   st->bot = newbot;
     933      268399 :   st->size = size;
     934      268399 :   return 1;
     935             : }
     936             : 
     937             : #else
     938             : #define PARI_STACK_ALIGN (0x40UL)
     939             : static void *
     940             : pari_mainstack_malloc(size_t s)
     941             : {
     942             :   return malloc(s); /* NOT pari_malloc, e_MEM would be deadly */
     943             : }
     944             : 
     945             : static void
     946             : pari_mainstack_mfree(void *s, size_t size) { (void) size; free(s); }
     947             : 
     948             : static int
     949             : pari_mainstack_setsize(struct pari_mainstack *st, size_t size)
     950             : {
     951             :   st->bot = st->top - size;
     952             :   st->size = size;
     953             :   return 1;
     954             : }
     955             : 
     956             : #endif
     957             : 
     958             : static const size_t MIN_STACK = 500032UL;
     959             : static size_t
     960      356652 : fix_size(size_t a)
     961             : {
     962      356652 :   size_t ps = PARI_STACK_ALIGN;
     963      356652 :   size_t b = a & ~(ps - 1); /* Align */
     964      356652 :   if (b < a && b < ~(ps - 1)) b += ps;
     965      356652 :   if (b < MIN_STACK) b = MIN_STACK;
     966      356652 :   return b;
     967             : }
     968             : 
     969             : static void
     970      178331 : pari_mainstack_alloc(int numerr, struct pari_mainstack *st, size_t rsize, size_t vsize)
     971             : {
     972      178331 :   size_t sizemax = vsize ? vsize: rsize, s = fix_size(sizemax);
     973             :   for (;;)
     974             :   {
     975      178331 :     st->vbot = (pari_sp)pari_mainstack_malloc(s);
     976      178331 :     if (st->vbot) break;
     977           0 :     if (s == MIN_STACK) pari_err(e_MEM); /* no way out. Die */
     978           0 :     s = fix_size(s >> 1);
     979           0 :     pari_warn(numerr, s);
     980             :   }
     981      178331 :   st->vsize = vsize ? s: 0;
     982      178331 :   st->rsize = minuu(rsize, s);
     983      178331 :   st->top = st->vbot+s;
     984      178331 :   if (!pari_mainstack_setsize(st, st->rsize))
     985             :   {
     986             :     /* This should never happen since we only decrease the allocated space */
     987           0 :     pari_err(e_MEM);
     988             :   }
     989      178331 :   st->memused = 0;
     990      178331 : }
     991             : 
     992             : static void
     993      178321 : pari_mainstack_free(struct pari_mainstack *st)
     994             : {
     995      178321 :   pari_mainstack_mfree((void*)st->vbot, st->vsize ? st->vsize : fix_size(st->rsize));
     996      178321 :   st->top = st->bot = st->vbot = 0;
     997      178321 :   st->size = st->vsize = 0;
     998      178321 : }
     999             : 
    1000             : static void
    1001         334 : pari_mainstack_resize(struct pari_mainstack *st, size_t rsize, size_t vsize)
    1002             : {
    1003         334 :   BLOCK_SIGINT_START;
    1004         334 :   pari_mainstack_free(st);
    1005         334 :   pari_mainstack_alloc(warnstack, st, rsize, vsize);
    1006         334 :   BLOCK_SIGINT_END;
    1007         334 : }
    1008             : 
    1009             : static void
    1010      177853 : pari_mainstack_use(struct pari_mainstack *st)
    1011             : {
    1012      177853 :   pari_mainstack = st;
    1013      177853 :   avma = st->top; /* don't use set_avma */
    1014      177853 : }
    1015             : 
    1016             : static void
    1017        1690 : paristack_alloc(size_t rsize, size_t vsize)
    1018             : {
    1019        1690 :   pari_mainstack_alloc(warnstack, pari_mainstack, rsize, vsize);
    1020        1690 :   pari_mainstack_use(pari_mainstack);
    1021        1690 : }
    1022             : 
    1023             : void
    1024           0 : paristack_setsize(size_t rsize, size_t vsize)
    1025             : {
    1026           0 :   pari_mainstack_resize(pari_mainstack, rsize, vsize);
    1027           0 :   pari_mainstack_use(pari_mainstack);
    1028           0 : }
    1029             : 
    1030             : void
    1031           0 : parivstack_resize(ulong newsize)
    1032             : {
    1033             :   size_t s;
    1034           0 :   if (newsize && newsize < pari_mainstack->rsize)
    1035           0 :     pari_err_DIM("stack sizes [parisizemax < parisize]");
    1036           0 :   if (newsize == pari_mainstack->vsize) return;
    1037           0 :   evalstate_reset();
    1038           0 :   paristack_setsize(pari_mainstack->rsize, newsize);
    1039           0 :   s = pari_mainstack->vsize ? pari_mainstack->vsize : pari_mainstack->rsize;
    1040           0 :   if (DEBUGMEM)
    1041           0 :     pari_warn(warner,"new maximum stack size = %lu (%.3f Mbytes)",
    1042             :               s, s/1048576.);
    1043           0 :   pari_init_errcatch();
    1044           0 :   cb_pari_err_recover(-1);
    1045             : }
    1046             : 
    1047             : void
    1048         340 : paristack_newrsize(ulong newsize)
    1049             : {
    1050         340 :   size_t s, vsize = pari_mainstack->vsize;
    1051         340 :   if (!newsize) newsize = pari_mainstack->rsize << 1;
    1052         340 :   if (newsize != pari_mainstack->rsize)
    1053         334 :     pari_mainstack_resize(pari_mainstack, newsize, vsize);
    1054         340 :   evalstate_reset();
    1055         340 :   s = pari_mainstack->rsize;
    1056         340 :   if (DEBUGMEM)
    1057         340 :     pari_warn(warner,"new stack size = %lu (%.3f Mbytes)", s, s/1048576.);
    1058         340 :   pari_init_errcatch();
    1059         340 :   cb_pari_err_recover(-1);
    1060           0 : }
    1061             : 
    1062             : void
    1063           0 : paristack_resize(ulong newsize)
    1064             : {
    1065           0 :   long size = pari_mainstack->size;
    1066           0 :   if (!newsize)
    1067           0 :     newsize = 2 * size;
    1068           0 :   newsize = minuu(newsize, pari_mainstack->vsize);
    1069           0 :   if (newsize <= pari_mainstack->size) return;
    1070           0 :   if (pari_mainstack_setsize(pari_mainstack, newsize))
    1071             :   {
    1072           0 :     if (DEBUGMEM)
    1073           0 :       pari_warn(warner, "increasing stack size to %lu", pari_mainstack->size);
    1074             :   }
    1075             :   else
    1076             :   {
    1077           0 :     pari_mainstack_setsize(pari_mainstack, size);
    1078           0 :     pari_err(e_STACK);
    1079             :   }
    1080             : }
    1081             : 
    1082             : void
    1083       90068 : parivstack_reset(void)
    1084             : {
    1085       90068 :   pari_mainstack_setsize(pari_mainstack, pari_mainstack->rsize);
    1086       90068 :   if (avma < pari_mainstack->bot)
    1087           0 :     pari_err_BUG("parivstack_reset [avma < bot]");
    1088       90068 : }
    1089             : 
    1090             : /* Enlarge the stack if needed such that the unused portion of the stack
    1091             :  * (between bot and avma) is large enough to contain x longs. */
    1092             : void
    1093           7 : new_chunk_resize(size_t x)
    1094             : {
    1095           7 :   if (pari_mainstack->vsize==0
    1096           7 :     || x > (avma-pari_mainstack->vbot) / sizeof(long)) pari_err(e_STACK);
    1097           0 :   while (x > (avma-pari_mainstack->bot) / sizeof(long))
    1098           0 :     paristack_resize(0);
    1099           0 : }
    1100             : 
    1101             : /*********************************************************************/
    1102             : /*                       PARI THREAD                                 */
    1103             : /*********************************************************************/
    1104             : 
    1105             : /* Initial PARI thread structure t with a stack of size s and
    1106             :  * argument arg */
    1107             : 
    1108             : static void
    1109      173581 : pari_thread_set_global(struct pari_global_state *gs)
    1110             : {
    1111      173581 :   push_localbitprec(gs->bitprec);
    1112      174943 :   pari_set_primetab(gs->primetab);
    1113      174904 :   pari_set_seadata(gs->seadata);
    1114      174857 :   pari_set_varstate(gs->varpriority, &gs->varstate);
    1115      175402 : }
    1116             : 
    1117             : static void
    1118      176307 : pari_thread_get_global(struct pari_global_state *gs)
    1119             : {
    1120      176307 :   gs->bitprec = get_localbitprec();
    1121      176307 :   gs->primetab = primetab;
    1122      176307 :   gs->seadata = pari_get_seadata();
    1123      176307 :   varstate_save(&gs->varstate);
    1124      176307 :   gs->varpriority = varpriority;
    1125      176307 : }
    1126             : 
    1127             : void
    1128      176307 : pari_thread_alloc(struct pari_thread *t, size_t s, GEN arg)
    1129             : {
    1130      176307 :   pari_mainstack_alloc(warnstackthread, &t->st,s,0);
    1131      176307 :   pari_thread_get_global(&t->gs);
    1132      176307 :   t->data = arg;
    1133      176307 : }
    1134             : 
    1135             : /* Initial PARI thread structure t with a stack of size s and virtual size v
    1136             :  * and argument arg */
    1137             : 
    1138             : void
    1139           0 : pari_thread_valloc(struct pari_thread *t, size_t s, size_t v, GEN arg)
    1140             : {
    1141           0 :   pari_mainstack_alloc(warnstackthread, &t->st,s,v);
    1142           0 :   pari_thread_get_global(&t->gs);
    1143           0 :   t->data = arg;
    1144           0 : }
    1145             : 
    1146             : void
    1147      176307 : pari_thread_free(struct pari_thread *t)
    1148             : {
    1149      176307 :   pari_mainstack_free(&t->st);
    1150      176307 : }
    1151             : 
    1152             : void
    1153      177807 : pari_thread_init(void)
    1154             : {
    1155             :   long var;
    1156      177807 :   pari_stackcheck_init((void*)&var);
    1157      177826 :   pari_init_blocks();
    1158      177746 :   pari_init_errcatch();
    1159      177691 :   pari_init_rand();
    1160      174581 :   pari_init_floats();
    1161      174929 :   pari_init_parser();
    1162      175453 :   pari_init_compiler();
    1163      175455 :   pari_init_evaluator();
    1164      175105 :   pari_init_files();
    1165      175313 : }
    1166             : 
    1167             : void
    1168      177118 : pari_thread_close(void)
    1169             : {
    1170      177118 :   pari_thread_close_files();
    1171      175693 :   pari_close_evaluator();
    1172      176843 :   pari_close_compiler();
    1173      173868 :   pari_close_parser();
    1174      176837 :   pari_close_floats();
    1175      174746 :   pari_close_blocks();
    1176      176448 : }
    1177             : 
    1178             : GEN
    1179      176181 : pari_thread_start(struct pari_thread *t)
    1180             : {
    1181      176181 :   pari_mainstack_use(&t->st);
    1182      176157 :   pari_thread_init();
    1183      173612 :   pari_thread_set_global(&t->gs);
    1184      175540 :   return t->data;
    1185             : }
    1186             : 
    1187             : /*********************************************************************/
    1188             : /*                       LIBPARI INIT / CLOSE                        */
    1189             : /*********************************************************************/
    1190             : 
    1191             : static void
    1192           0 : pari_exit(void)
    1193             : {
    1194           0 :   err_printf("  ***   Error in the PARI system. End of program.\n");
    1195           0 :   exit(1);
    1196             : }
    1197             : 
    1198             : static void
    1199           0 : dflt_err_recover(long errnum) { (void) errnum; pari_exit(); }
    1200             : 
    1201             : static void
    1202           0 : dflt_pari_quit(long err) { (void)err; /*do nothing*/; }
    1203             : 
    1204             : static int pari_err_display(GEN err);
    1205             : 
    1206             : /* initialize PARI data. Initialize [new|old]fun to NULL for default set. */
    1207             : void
    1208        1690 : pari_init_opts(size_t parisize, ulong maxprime, ulong init_opts)
    1209             : {
    1210             :   ulong u;
    1211             : 
    1212        1690 :   pari_mt_nbthreads = 0;
    1213        1690 :   cb_pari_quit = dflt_pari_quit;
    1214        1690 :   cb_pari_init_histfile = NULL;
    1215        1690 :   cb_pari_get_line_interactive = NULL;
    1216        1690 :   cb_pari_fgets_interactive = NULL;
    1217        1690 :   cb_pari_whatnow = NULL;
    1218        1690 :   cb_pari_handle_exception = NULL;
    1219        1690 :   cb_pari_err_handle = pari_err_display;
    1220        1690 :   cb_pari_pre_recover = NULL;
    1221        1690 :   cb_pari_break_loop = NULL;
    1222        1690 :   cb_pari_is_interactive = NULL;
    1223        1690 :   cb_pari_start_output = NULL;
    1224        1690 :   cb_pari_sigint = dflt_sigint_fun;
    1225        1690 :   if (init_opts&INIT_JMPm) cb_pari_err_recover = dflt_err_recover;
    1226             : 
    1227        1690 :   pari_stackcheck_init(&u);
    1228        1690 :   pari_init_homedir();
    1229        1690 :   if (init_opts&INIT_DFTm) {
    1230           0 :     pari_init_defaults();
    1231           0 :     GP_DATA = default_gp_data();
    1232           0 :     pari_init_paths();
    1233             :   }
    1234             : 
    1235        1690 :   pari_mainstack = (struct pari_mainstack *) malloc(sizeof(*pari_mainstack));
    1236        1690 :   paristack_alloc(parisize, 0);
    1237        1690 :   init_universal_constants();
    1238        1690 :   diffptr = NULL;
    1239        1690 :   if (!(init_opts&INIT_noPRIMEm))
    1240             :   {
    1241           0 :     GP_DATA->primelimit = maxprime;
    1242           0 :     pari_init_primes(GP_DATA->primelimit);
    1243             :   }
    1244        1690 :   if (!(init_opts&INIT_noINTGMPm)) pari_kernel_init();
    1245        1690 :   pari_init_graphics();
    1246        1690 :   pari_thread_init();
    1247        1690 :   pari_set_primetab(NULL);
    1248        1690 :   pari_set_seadata(NULL);
    1249        1690 :   pari_init_functions();
    1250        1690 :   pari_init_export();
    1251        1690 :   pari_var_init();
    1252        1690 :   pari_init_timer();
    1253        1690 :   pari_init_buffers();
    1254        1690 :   (void)getabstime();
    1255        1690 :   try_to_recover = 1;
    1256        1690 :   if (!(init_opts&INIT_noIMTm)) pari_mt_init();
    1257        1690 :   if ((init_opts&INIT_SIGm)) pari_sig_init(pari_sighandler);
    1258        1690 : }
    1259             : 
    1260             : void
    1261           0 : pari_init(size_t parisize, ulong maxprime)
    1262           0 : { pari_init_opts(parisize, maxprime, INIT_JMPm | INIT_SIGm | INIT_DFTm); }
    1263             : 
    1264             : void
    1265        1680 : pari_close_opts(ulong init_opts)
    1266             : {
    1267             :   long i;
    1268             : 
    1269        1680 :   BLOCK_SIGINT_START;
    1270        1680 :   if ((init_opts&INIT_SIGm)) pari_sig_init(SIG_DFL);
    1271        1680 :   if (!(init_opts&INIT_noIMTm)) pari_mt_close();
    1272             : 
    1273        1680 :   pari_var_close(); /* must come before destruction of functions_hash */
    1274      228480 :   for (i = 0; i < functions_tblsz; i++)
    1275             :   {
    1276      226800 :     entree *ep = functions_hash[i];
    1277     2336331 :     while (ep) {
    1278     2109531 :       entree *EP = ep->next;
    1279     2109531 :       if (!EpSTATIC(ep)) { freeep(ep); free(ep); }
    1280     2109531 :       ep = EP;
    1281             :     }
    1282             :   }
    1283        1680 :   pari_close_mf();
    1284        1680 :   pari_thread_close();
    1285        1680 :   pari_close_export();
    1286        1680 :   pari_close_files();
    1287        1680 :   pari_close_homedir();
    1288        1680 :   if (!(init_opts&INIT_noINTGMPm)) pari_kernel_close();
    1289             : 
    1290        1680 :   free((void*)functions_hash);
    1291        1680 :   free((void*)defaults_hash);
    1292        1680 :   if (diffptr) pari_close_primes();
    1293        1680 :   free(current_logfile);
    1294        1680 :   free(current_psfile);
    1295        1680 :   pari_mainstack_free(pari_mainstack);
    1296        1680 :   free((void*)pari_mainstack);
    1297        1680 :   pari_stack_delete(&s_MODULES);
    1298        1680 :   if (pari_datadir) free(pari_datadir);
    1299        1680 :   if (init_opts&INIT_DFTm)
    1300             :   { /* delete GP_DATA */
    1301        1680 :     pari_close_paths();
    1302        1680 :     if (GP_DATA->hist->v) free((void*)GP_DATA->hist->v);
    1303        1680 :     if (GP_DATA->pp->cmd) free((void*)GP_DATA->pp->cmd);
    1304        1680 :     if (GP_DATA->help) free((void*)GP_DATA->help);
    1305        1680 :     if (GP_DATA->plothsizes) free((void*)GP_DATA->plothsizes);
    1306        1680 :     if (GP_DATA->colormap) pari_free(GP_DATA->colormap);
    1307        1680 :     if (GP_DATA->graphcolors) pari_free(GP_DATA->graphcolors);
    1308        1680 :     free((void*)GP_DATA->prompt);
    1309        1680 :     free((void*)GP_DATA->prompt_cont);
    1310        1680 :     free((void*)GP_DATA->histfile);
    1311             :   }
    1312        1680 :   BLOCK_SIGINT_END;
    1313        1680 : }
    1314             : 
    1315             : void
    1316        1680 : pari_close(void)
    1317        1680 : { pari_close_opts(INIT_JMPm | INIT_SIGm | INIT_DFTm); }
    1318             : 
    1319             : /*******************************************************************/
    1320             : /*                                                                 */
    1321             : /*                         ERROR RECOVERY                          */
    1322             : /*                                                                 */
    1323             : /*******************************************************************/
    1324             : void
    1325      108908 : gp_context_save(struct gp_context* rec)
    1326             : {
    1327      108908 :   rec->prettyp = GP_DATA->fmt->prettyp;
    1328      108908 :   rec->listloc = next_block;
    1329      108908 :   rec->iferr_env = iferr_env;
    1330      108908 :   rec->err_data  = global_err_data;
    1331      108908 :   varstate_save(&rec->var);
    1332      108908 :   evalstate_save(&rec->eval);
    1333      108908 :   parsestate_save(&rec->parse);
    1334      108908 :   filestate_save(&rec->file);
    1335      108908 : }
    1336             : 
    1337             : void
    1338        9974 : gp_context_restore(struct gp_context* rec)
    1339             : {
    1340             :   long i;
    1341             : 
    1342        9974 :   if (!try_to_recover) return;
    1343             :   /* disable gp_context_restore() and SIGINT */
    1344        9974 :   try_to_recover = 0;
    1345        9974 :   BLOCK_SIGINT_START
    1346        9974 :   if (DEBUGMEM>2) err_printf("entering recover(), loc = %ld\n", rec->listloc);
    1347        9974 :   evalstate_restore(&rec->eval);
    1348        9974 :   parsestate_restore(&rec->parse);
    1349        9974 :   filestate_restore(&rec->file);
    1350        9974 :   global_err_data = rec->err_data;
    1351        9974 :   iferr_env = rec->iferr_env;
    1352        9974 :   GP_DATA->fmt->prettyp = rec->prettyp;
    1353             : 
    1354     1356464 :   for (i = 0; i < functions_tblsz; i++)
    1355             :   {
    1356     1346490 :     entree *ep = functions_hash[i];
    1357    14543005 :     while (ep)
    1358             :     {
    1359    13196515 :       entree *EP = ep->next;
    1360    13196515 :       switch(EpVALENCE(ep))
    1361             :       {
    1362      219334 :         case EpVAR:
    1363      219684 :           while (pop_val_if_newer(ep,rec->listloc)) /* empty */;
    1364      219334 :           break;
    1365      629357 :         case EpNEW: break;
    1366             :       }
    1367    13196515 :       ep = EP;
    1368             :     }
    1369             :   }
    1370        9974 :   varstate_restore(&rec->var);
    1371        9974 :   if (DEBUGMEM>2) err_printf("leaving recover()\n");
    1372        9974 :   BLOCK_SIGINT_END
    1373        9974 :   try_to_recover = 1;
    1374             : }
    1375             : 
    1376             : static void
    1377        9899 : err_recover(long numerr)
    1378             : {
    1379        9899 :   if (cb_pari_pre_recover)
    1380        9899 :     cb_pari_pre_recover(numerr);
    1381           0 :   evalstate_reset();
    1382           0 :   killallfiles();
    1383           0 :   pari_init_errcatch();
    1384           0 :   cb_pari_err_recover(numerr);
    1385           0 : }
    1386             : 
    1387             : static void
    1388       10526 : err_init(void)
    1389             : {
    1390             :   /* make sure pari_err msg starts at the beginning of line */
    1391       10526 :   if (!pari_last_was_newline()) pari_putc('\n');
    1392       10526 :   pariOut->flush();
    1393       10526 :   pariErr->flush();
    1394       10526 :   out_term_color(pariErr, c_ERR);
    1395       10526 : }
    1396             : 
    1397             : static void
    1398       10421 : err_init_msg(int user)
    1399             : {
    1400             :   const char *gp_function_name;
    1401       10421 :   out_puts(pariErr, "  *** ");
    1402       10421 :   if (!user && (gp_function_name = closure_func_err()))
    1403        6883 :     out_printf(pariErr, "%s: ", gp_function_name);
    1404             :   else
    1405        3538 :     out_puts(pariErr, "  ");
    1406       10421 : }
    1407             : 
    1408             : void
    1409         606 : pari_warn(int numerr, ...)
    1410             : {
    1411             :   char *ch1;
    1412             :   va_list ap;
    1413             : 
    1414         606 :   va_start(ap,numerr);
    1415             : 
    1416         606 :   err_init();
    1417         606 :   err_init_msg(numerr==warnuser || numerr==warnstack);
    1418         606 :   switch (numerr)
    1419             :   {
    1420           7 :     case warnuser:
    1421           7 :       out_puts(pariErr, "user warning: ");
    1422           7 :       out_print0(pariErr, NULL, va_arg(ap, GEN), f_RAW);
    1423           7 :       break;
    1424             : 
    1425           0 :     case warnmem:
    1426           0 :       out_puts(pariErr, "collecting garbage in "); ch1=va_arg(ap, char*);
    1427           0 :       out_vprintf(pariErr, ch1,ap); out_putc(pariErr, '.');
    1428           0 :       break;
    1429             : 
    1430         599 :     case warner:
    1431         599 :       out_puts(pariErr, "Warning: "); ch1=va_arg(ap, char*);
    1432         599 :       out_vprintf(pariErr, ch1,ap); out_putc(pariErr, '.');
    1433         599 :       break;
    1434             : 
    1435           0 :     case warnprec:
    1436           0 :       out_vprintf(pariErr, "Warning: increasing prec in %s; new prec = %ld",
    1437             :                       ap);
    1438           0 :       break;
    1439             : 
    1440           0 :     case warnfile:
    1441           0 :       out_puts(pariErr, "Warning: failed to "),
    1442           0 :       ch1 = va_arg(ap, char*);
    1443           0 :       out_printf(pariErr, "%s: %s", ch1, va_arg(ap, char*));
    1444           0 :       break;
    1445             : 
    1446           0 :     case warnstack:
    1447             :     case warnstackthread:
    1448             :     {
    1449           0 :       ulong  s = va_arg(ap, ulong);
    1450             :       char buf[128];
    1451           0 :       const char * stk = numerr == warnstackthread
    1452           0 :                          || mt_is_thread() ? "thread": "PARI";
    1453           0 :       sprintf(buf,"Warning: not enough memory, new %s stack %lu", stk, s);
    1454           0 :       out_puts(pariErr,buf);
    1455           0 :       break;
    1456             :     }
    1457             :   }
    1458         606 :   va_end(ap);
    1459         606 :   out_term_color(pariErr, c_NONE);
    1460         606 :   out_putc(pariErr, '\n');
    1461         606 :   pariErr->flush();
    1462         606 : }
    1463             : void
    1464           0 : pari_sigint(const char *time_s)
    1465             : {
    1466           0 :   int recover=0;
    1467           0 :   BLOCK_SIGALRM_START
    1468           0 :   err_init();
    1469           0 :   closure_err(0);
    1470           0 :   err_init_msg(0);
    1471           0 :   out_puts(pariErr, "user interrupt after ");
    1472           0 :   out_puts(pariErr, time_s);
    1473           0 :   out_term_color(pariErr, c_NONE);
    1474           0 :   pariErr->flush();
    1475           0 :   if (cb_pari_handle_exception)
    1476           0 :     recover = cb_pari_handle_exception(-1);
    1477           0 :   if (!recover && !block)
    1478           0 :     PARI_SIGINT_pending = 0;
    1479           0 :   BLOCK_SIGINT_END
    1480           0 :   if (!recover) err_recover(e_MISC);
    1481           0 : }
    1482             : 
    1483             : #define retmkerr2(x,y)\
    1484             :   do { GEN _v = cgetg(3, t_ERROR);\
    1485             :        _v[1] = (x);\
    1486             :        gel(_v,2) = (y); return _v; } while(0)
    1487             : #define retmkerr3(x,y,z)\
    1488             :   do { GEN _v = cgetg(4, t_ERROR);\
    1489             :        _v[1] = (x);\
    1490             :        gel(_v,2) = (y);\
    1491             :        gel(_v,3) = (z); return _v; } while(0)
    1492             : #define retmkerr4(x,y,z,t)\
    1493             :   do { GEN _v = cgetg(5, t_ERROR);\
    1494             :        _v[1] = (x);\
    1495             :        gel(_v,2) = (y);\
    1496             :        gel(_v,3) = (z);\
    1497             :        gel(_v,4) = (t); return _v; } while(0)
    1498             : #define retmkerr5(x,y,z,t,u)\
    1499             :   do { GEN _v = cgetg(6, t_ERROR);\
    1500             :        _v[1] = (x);\
    1501             :        gel(_v,2) = (y);\
    1502             :        gel(_v,3) = (z);\
    1503             :        gel(_v,4) = (t);\
    1504             :        gel(_v,5) = (u); return _v; } while(0)
    1505             : #define retmkerr6(x,y,z,t,u,v)\
    1506             :   do { GEN _v = cgetg(7, t_ERROR);\
    1507             :        _v[1] = (x);\
    1508             :        gel(_v,2) = (y);\
    1509             :        gel(_v,3) = (z);\
    1510             :        gel(_v,4) = (t);\
    1511             :        gel(_v,5) = (u);\
    1512             :        gel(_v,6) = (v); return _v; } while(0)
    1513             : 
    1514             : static GEN
    1515       45012 : pari_err2GEN(long numerr, va_list ap)
    1516             : {
    1517       45012 :   switch ((enum err_list) numerr)
    1518             :   {
    1519         105 :   case e_SYNTAX:
    1520             :     {
    1521         105 :       const char *msg = va_arg(ap, char*);
    1522         105 :       const char *s = va_arg(ap,char *);
    1523         105 :       const char *entry = va_arg(ap,char *);
    1524         105 :       retmkerr3(numerr,strtoGENstr(msg), mkvecsmall2((long)s,(long)entry));
    1525             :     }
    1526         191 :   case e_MISC: case e_ALARM:
    1527             :     {
    1528         191 :       const char *ch1 = va_arg(ap, char*);
    1529         191 :       retmkerr2(numerr, gvsprintf(ch1,ap));
    1530             :     }
    1531        2737 :   case e_NOTFUNC:
    1532             :   case e_USER:
    1533        2737 :     retmkerr2(numerr,va_arg(ap, GEN));
    1534           0 :   case e_FILE:
    1535             :     {
    1536           0 :       const char *f = va_arg(ap, const char*);
    1537           0 :       retmkerr3(numerr, strtoGENstr(f), strtoGENstr(va_arg(ap, char*)));
    1538             :     }
    1539          36 :   case e_FILEDESC:
    1540             :     {
    1541          36 :       const char *f = va_arg(ap, const char*);
    1542          36 :       retmkerr3(numerr, strtoGENstr(f), stoi(va_arg(ap, long)));
    1543             :     }
    1544        1484 :   case e_OVERFLOW:
    1545             :   case e_IMPL:
    1546             :   case e_DIM:
    1547             :   case e_CONSTPOL:
    1548             :   case e_ROOTS0:
    1549             :   case e_FLAG:
    1550             :   case e_PREC:
    1551             :   case e_BUG:
    1552             :   case e_ARCH:
    1553             :   case e_PACKAGE:
    1554        1484 :     retmkerr2(numerr, strtoGENstr(va_arg(ap, char*)));
    1555        1029 :   case e_MODULUS:
    1556             :   case e_VAR:
    1557             :     {
    1558        1029 :       const char *f = va_arg(ap, const char*);
    1559        1029 :       GEN x = va_arg(ap, GEN);
    1560        1029 :       GEN y = va_arg(ap, GEN);
    1561        1029 :       retmkerr4(numerr, strtoGENstr(f), x,y);
    1562             :     }
    1563       33043 :   case e_INV:
    1564             :   case e_IRREDPOL:
    1565             :   case e_PRIME:
    1566             :   case e_SQRTN:
    1567             :   case e_TYPE:
    1568             :     {
    1569       33043 :       const char *f = va_arg(ap, const char*);
    1570       33043 :       GEN x = va_arg(ap, GEN);
    1571       33043 :       retmkerr3(numerr, strtoGENstr(f), x);
    1572             :     }
    1573        3822 :   case e_COPRIME: case e_OP: case e_TYPE2:
    1574             :     {
    1575        3822 :       const char *f = va_arg(ap, const char*);
    1576        3822 :       GEN x = va_arg(ap, GEN);
    1577        3822 :       GEN y = va_arg(ap, GEN);
    1578        3822 :       retmkerr4(numerr,strtoGENstr(f),x,y);
    1579             :     }
    1580         200 :   case e_COMPONENT:
    1581             :     {
    1582         200 :       const char *f= va_arg(ap, const char *);
    1583         200 :       const char *op = va_arg(ap, const char *);
    1584         200 :       GEN l = va_arg(ap, GEN);
    1585         200 :       GEN x = va_arg(ap, GEN);
    1586         200 :       retmkerr5(numerr,strtoGENstr(f),strtoGENstr(op),l,x);
    1587             :     }
    1588        2162 :   case e_DOMAIN:
    1589             :     {
    1590        2162 :       const char *f = va_arg(ap, const char*);
    1591        2162 :       const char *v = va_arg(ap, const char *);
    1592        2162 :       const char *op = va_arg(ap, const char *);
    1593        2162 :       GEN l = va_arg(ap, GEN);
    1594        2162 :       GEN x = va_arg(ap, GEN);
    1595        2162 :       retmkerr6(numerr,strtoGENstr(f),strtoGENstr(v),strtoGENstr(op),l,x);
    1596             :     }
    1597         196 :   case e_PRIORITY:
    1598             :     {
    1599         196 :       const char *f = va_arg(ap, const char*);
    1600         196 :       GEN x = va_arg(ap, GEN);
    1601         196 :       const char *op = va_arg(ap, const char *);
    1602         196 :       long v = va_arg(ap, long);
    1603         196 :       retmkerr5(numerr,strtoGENstr(f),x,strtoGENstr(op),stoi(v));
    1604             :     }
    1605           0 :   case e_MAXPRIME:
    1606           0 :     retmkerr2(numerr, utoi(va_arg(ap, ulong)));
    1607           7 :   case e_STACK:
    1608           7 :     return err_e_STACK;
    1609           0 :   case e_STACKTHREAD:
    1610           0 :     retmkerr3(numerr, utoi(va_arg(ap, ulong)), utoi(va_arg(ap, ulong)));
    1611           0 :   default:
    1612           0 :     return mkerr(numerr);
    1613             :   }
    1614             : }
    1615             : 
    1616             : static char *
    1617        7028 : type_dim(GEN x)
    1618             : {
    1619        7028 :   char *v = stack_malloc(64);
    1620        7028 :   switch(typ(x))
    1621             :   {
    1622          84 :     case t_MAT:
    1623             :     {
    1624          84 :       long l = lg(x), r = (l == 1)? 1: lgcols(x);
    1625          84 :       sprintf(v, "t_MAT (%ldx%ld)", r-1,l-1);
    1626          84 :       break;
    1627             :     }
    1628          91 :     case t_COL:
    1629          91 :       sprintf(v, "t_COL (%ld elts)", lg(x)-1);
    1630          91 :       break;
    1631         133 :     case t_VEC:
    1632         133 :       sprintf(v, "t_VEC (%ld elts)", lg(x)-1);
    1633         133 :       break;
    1634        6720 :     default:
    1635        6720 :       v = (char*)type_name(typ(x));
    1636             :   }
    1637        7028 :   return v;
    1638             : }
    1639             : 
    1640             : static char *
    1641        2443 : gdisplay(GEN x)
    1642             : {
    1643        2443 :   char *s = GENtostr_raw(x);
    1644        2443 :   if (strlen(s) < 1600) return s;
    1645          21 :   if (! GP_DATA->breakloop) return (char*)"(...)";
    1646           0 :   return stack_sprintf("\n  ***  (...) Huge %s omitted; you can access it via dbg_err()", type_name(typ(x)));
    1647             : }
    1648             : 
    1649             : char *
    1650       17515 : pari_err2str(GEN e)
    1651             : {
    1652       17515 :   long numerr = err_get_num(e);
    1653       17515 :   switch ((enum err_list) numerr)
    1654             :   {
    1655           0 :   case e_ALARM:
    1656           0 :     return pari_sprintf("alarm interrupt after %Ps.",gel(e,2));
    1657         182 :   case e_MISC:
    1658         182 :     return pari_sprintf("%Ps.",gel(e,2));
    1659             : 
    1660           0 :   case e_ARCH:
    1661           0 :     return pari_sprintf("sorry, '%Ps' not available on this system.",gel(e,2));
    1662          14 :   case e_BUG:
    1663          14 :     return pari_sprintf("bug in %Ps, please report.",gel(e,2));
    1664           7 :   case e_CONSTPOL:
    1665           7 :     return pari_sprintf("constant polynomial in %Ps.", gel(e,2));
    1666          70 :   case e_COPRIME:
    1667         210 :     return pari_sprintf("elements not coprime in %Ps:\n    %s\n    %s",
    1668          70 :                         gel(e,2), gdisplay(gel(e,3)), gdisplay(gel(e,4)));
    1669         529 :   case e_DIM:
    1670         529 :     return pari_sprintf("inconsistent dimensions in %Ps.", gel(e,2));
    1671           0 :   case e_FILE:
    1672           0 :     return pari_sprintf("error opening %Ps: `%Ps'.", gel(e,2), gel(e,3));
    1673          36 :   case e_FILEDESC:
    1674          36 :     return pari_sprintf("invalid file descriptor in %Ps [%Ps]", gel(e,2), gel(e,3));
    1675          84 :   case e_FLAG:
    1676          84 :     return pari_sprintf("invalid flag in %Ps.", gel(e,2));
    1677         420 :   case e_IMPL:
    1678         420 :     return pari_sprintf("sorry, %Ps is not yet implemented.", gel(e,2));
    1679           0 :   case e_PACKAGE:
    1680           0 :     return pari_sprintf("package %Ps is required, please install it.", gel(e,2));
    1681         595 :   case e_INV:
    1682         595 :     return pari_sprintf("impossible inverse in %Ps: %s.", gel(e,2),
    1683         595 :                         gdisplay(gel(e,3)));
    1684          49 :   case e_IRREDPOL:
    1685          98 :     return pari_sprintf("not an irreducible polynomial in %Ps: %s.",
    1686          49 :                         gel(e,2), gdisplay(gel(e,3)));
    1687           0 :   case e_MAXPRIME:
    1688             :     {
    1689           0 :       const char * msg = "not enough precomputed primes";
    1690           0 :       ulong c = itou(gel(e,2));
    1691           0 :       if (c) return pari_sprintf("%s, need primelimit ~ %lu.",msg, c);
    1692           0 :       else   return pari_strdup(msg);
    1693             :     }
    1694           0 :   case e_MEM:
    1695           0 :     return pari_strdup("not enough memory");
    1696         749 :   case e_MODULUS:
    1697             :     {
    1698         749 :       GEN x = gel(e,3), y = gel(e,4);
    1699         749 :       return pari_sprintf("inconsistent moduli in %Ps: %s != %s",
    1700         749 :                           gel(e,2), gdisplay(x), gdisplay(y));
    1701             :     }
    1702           0 :   case e_NONE: return NULL;
    1703        2723 :   case e_NOTFUNC:
    1704        2723 :     return pari_strdup("not a function in function call");
    1705        3514 :   case e_OP: case e_TYPE2:
    1706             :     {
    1707        3514 :       pari_sp av = avma;
    1708             :       char *v;
    1709        3514 :       const char *f, *op = GSTR(gel(e,2));
    1710        3514 :       const char *what = numerr == e_OP? "inconsistent": "forbidden";
    1711        3514 :       GEN x = gel(e,3);
    1712        3514 :       GEN y = gel(e,4);
    1713        3514 :       switch(*op)
    1714             :       {
    1715          14 :       case '+': f = "addition"; break;
    1716          91 :       case '*': f = "multiplication"; break;
    1717        2744 :       case '/': case '%': case '\\': f = "division"; break;
    1718           0 :       case '=': op = "-->"; f = "assignment"; break;
    1719         665 :       default:  f = op; op = ","; break;
    1720             :       }
    1721        3514 :       v = pari_sprintf("%s %s %s %s %s.", what,f,type_dim(x),op,type_dim(y));
    1722        3514 :       set_avma(av); return v;
    1723             :     }
    1724         200 :   case e_COMPONENT:
    1725             :     {
    1726         200 :       const char *f= GSTR(gel(e,2));
    1727         200 :       const char *op= GSTR(gel(e,3));
    1728         200 :       GEN l = gel(e,4);
    1729         200 :       if (!*f)
    1730         147 :         return pari_sprintf("non-existent component: index %s %Ps",op,l);
    1731          53 :       return pari_sprintf("non-existent component in %s: index %s %Ps",f,op,l);
    1732             :     }
    1733        2052 :   case e_DOMAIN:
    1734             :     {
    1735        2052 :       const char *f = GSTR(gel(e,2));
    1736        2052 :       const char *v = GSTR(gel(e,3));
    1737        2052 :       const char *op= GSTR(gel(e,4));
    1738        2052 :       GEN l = gel(e,5);
    1739        2052 :       if (!*op)
    1740          28 :         return pari_sprintf("domain error in %s: %s out of range",f,v);
    1741        2024 :       return pari_sprintf("domain error in %s: %s %s %Ps",f,v,op,l);
    1742             :     }
    1743         147 :   case e_PRIORITY:
    1744             :     {
    1745         147 :       const char *f = GSTR(gel(e,2));
    1746         147 :       long vx = gvar(gel(e,3));
    1747         147 :       const char *op= GSTR(gel(e,4));
    1748         147 :       long v = itos(gel(e,5));
    1749         147 :       return pari_sprintf("incorrect priority in %s: variable %Ps %s %Ps",f,
    1750             :              pol_x(vx), op, pol_x(v));
    1751             :     }
    1752         119 :   case e_OVERFLOW:
    1753         119 :     return pari_sprintf("overflow in %Ps.", gel(e,2));
    1754         217 :   case e_PREC:
    1755         217 :     return pari_sprintf("precision too low in %Ps.", gel(e,2));
    1756          77 :   case e_PRIME:
    1757         154 :     return pari_sprintf("not a prime number in %Ps: %s.",
    1758          77 :                         gel(e,2), gdisplay(gel(e,3)));
    1759          49 :   case e_ROOTS0:
    1760          49 :     return pari_sprintf("zero polynomial in %Ps.", gel(e,2));
    1761          84 :   case e_SQRTN:
    1762         168 :     return pari_sprintf("not an n-th power residue in %Ps: %s.",
    1763          84 :                         gel(e,2), gdisplay(gel(e,3)));
    1764           7 :   case e_STACK:
    1765             :   case e_STACKTHREAD:
    1766             :     {
    1767           7 :       const char *stack = numerr == e_STACK? "PARI": "thread";
    1768           7 :       const char *var = numerr == e_STACK? "parisizemax": "threadsizemax";
    1769           0 :       size_t rsize = numerr == e_STACKTHREAD && GP_DATA->threadsize ?
    1770           7 :                                 GP_DATA->threadsize: pari_mainstack->rsize;
    1771           7 :       size_t vsize = numerr == e_STACK? pari_mainstack->vsize:
    1772           0 :                                         GP_DATA->threadsizemax;
    1773           7 :       char *buf = (char *) pari_malloc(512*sizeof(char));
    1774           7 :       if (vsize)
    1775             :       {
    1776           0 :         sprintf(buf, "the %s stack overflows !\n"
    1777             :             "  current stack size: %lu (%.3f Mbytes)\n"
    1778             :             "  [hint] you can increase '%s' using default()\n",
    1779           0 :             stack, (ulong)vsize, (double)vsize/1048576., var);
    1780             :       }
    1781             :       else
    1782             :       {
    1783           7 :         sprintf(buf, "the %s stack overflows !\n"
    1784             :             "  current stack size: %lu (%.3f Mbytes)\n"
    1785             :             "  [hint] set '%s' to a non-zero value in your GPRC\n",
    1786           7 :             stack, (ulong)rsize, (double)rsize/1048576., var);
    1787             :       }
    1788           7 :       return buf;
    1789             :     }
    1790           0 :   case e_SYNTAX:
    1791           0 :     return pari_strdup(GSTR(gel(e,2)));
    1792        5297 :   case e_TYPE:
    1793       10594 :     return pari_sprintf("incorrect type in %Ps (%s).",
    1794        5297 :                         gel(e,2), type_name(typ(gel(e,3))));
    1795          14 :   case e_USER:
    1796          14 :     return pari_sprint0("user error: ", gel(e,2), f_RAW);
    1797         280 :   case e_VAR:
    1798             :     {
    1799         280 :       GEN x = gel(e,3), y = gel(e,4);
    1800         840 :       return pari_sprintf("inconsistent variables in %Ps, %Ps != %Ps.",
    1801         280 :                           gel(e,2), pol_x(varn(x)), pol_x(varn(y)));
    1802             :     }
    1803             :   }
    1804             :   return NULL; /*LCOV_EXCL_LINE*/
    1805             : }
    1806             : 
    1807             : static int
    1808        9920 : pari_err_display(GEN err)
    1809             : {
    1810        9920 :   long numerr=err_get_num(err);
    1811        9920 :   err_init();
    1812        9920 :   if (numerr==e_SYNTAX)
    1813             :   {
    1814         105 :     const char *msg = GSTR(gel(err,2));
    1815         105 :     const char *s     = (const char *) gmael(err,3,1);
    1816         105 :     const char *entry = (const char *) gmael(err,3,2);
    1817         105 :     print_errcontext(pariErr, msg, s, entry);
    1818             :   }
    1819             :   else
    1820             :   {
    1821             :     char *s;
    1822        9815 :     closure_err(0);
    1823        9815 :     err_init_msg(numerr==e_USER);
    1824        9815 :     s = pari_err2str(err); pariErr->puts(s); pari_free(s);
    1825        9815 :     if (numerr==e_NOTFUNC)
    1826             :     {
    1827        2723 :       GEN fun = gel(err,2);
    1828        2723 :       if (gequalX(fun))
    1829             :       {
    1830        2723 :         entree *ep = varentries[varn(fun)];
    1831        2723 :         const char *t = ep->name;
    1832        2723 :         if (cb_pari_whatnow) cb_pari_whatnow(pariErr,t,1);
    1833             :       }
    1834             :     }
    1835             :   }
    1836        9906 :   out_term_color(pariErr, c_NONE);
    1837        9906 :   pariErr->flush(); return 0;
    1838             : }
    1839             : 
    1840             : void
    1841       45031 : pari_err(int numerr, ...)
    1842             : {
    1843             :   va_list ap;
    1844             :   GEN E;
    1845             : 
    1846       45031 :   va_start(ap,numerr);
    1847             : 
    1848       45031 :   if (numerr)
    1849       45012 :     E = pari_err2GEN(numerr,ap);
    1850             :   else
    1851             :   {
    1852          19 :     E = va_arg(ap,GEN);
    1853          19 :     numerr = err_get_num(E);
    1854             :   }
    1855       45029 :   global_err_data = E;
    1856       45029 :   if (*iferr_env) longjmp(*iferr_env, numerr);
    1857        9926 :   mt_err_recover(numerr);
    1858        9920 :   va_end(ap);
    1859       19826 :   if (cb_pari_err_handle &&
    1860        9920 :       cb_pari_err_handle(E)) return;
    1861       19803 :   if (cb_pari_handle_exception &&
    1862        9904 :       cb_pari_handle_exception(numerr)) return;
    1863        9899 :   err_recover(numerr);
    1864             : }
    1865             : 
    1866             : GEN
    1867       70203 : pari_err_last(void) { return global_err_data; }
    1868             : 
    1869             : const char *
    1870       26927 : numerr_name(long numerr)
    1871             : {
    1872       26927 :   switch ((enum err_list) numerr)
    1873             :   {
    1874           0 :   case e_ALARM:    return "e_ALARM";
    1875           0 :   case e_ARCH:     return "e_ARCH";
    1876           0 :   case e_BUG:      return "e_BUG";
    1877           0 :   case e_COMPONENT: return "e_COMPONENT";
    1878           0 :   case e_CONSTPOL: return "e_CONSTPOL";
    1879           0 :   case e_COPRIME:  return "e_COPRIME";
    1880           0 :   case e_DIM:      return "e_DIM";
    1881          56 :   case e_DOMAIN:   return "e_DOMAIN";
    1882           0 :   case e_FILE:     return "e_FILE";
    1883           0 :   case e_FILEDESC: return "e_FILEDESC";
    1884           7 :   case e_FLAG:     return "e_FLAG";
    1885          28 :   case e_IMPL:     return "e_IMPL";
    1886       19101 :   case e_INV:      return "e_INV";
    1887           0 :   case e_IRREDPOL: return "e_IRREDPOL";
    1888           0 :   case e_MAXPRIME: return "e_MAXPRIME";
    1889           0 :   case e_MEM:      return "e_MEM";
    1890           0 :   case e_MISC:     return "e_MISC";
    1891           0 :   case e_MODULUS:  return "e_MODULUS";
    1892           0 :   case e_NONE:     return "e_NONE";
    1893           0 :   case e_NOTFUNC:  return "e_NOTFUNC";
    1894           0 :   case e_OP:       return "e_OP";
    1895           0 :   case e_OVERFLOW: return "e_OVERFLOW";
    1896           0 :   case e_PACKAGE:  return "e_PACKAGE";
    1897           0 :   case e_PREC:     return "e_PREC";
    1898           0 :   case e_PRIME:    return "e_PRIME";
    1899          49 :   case e_PRIORITY: return "e_PRIORITY";
    1900           0 :   case e_ROOTS0:   return "e_ROOTS0";
    1901           0 :   case e_SQRTN:    return "e_SQRTN";
    1902           0 :   case e_STACK:    return "e_STACK";
    1903           0 :   case e_SYNTAX:   return "e_SYNTAX";
    1904           0 :   case e_STACKTHREAD:   return "e_STACKTHREAD";
    1905           0 :   case e_TYPE2:    return "e_TYPE2";
    1906        7686 :   case e_TYPE:     return "e_TYPE";
    1907           0 :   case e_USER:     return "e_USER";
    1908           0 :   case e_VAR:      return "e_VAR";
    1909             :   }
    1910           0 :   return "invalid error number";
    1911             : }
    1912             : 
    1913             : long
    1914          21 : name_numerr(const char *s)
    1915             : {
    1916          21 :   if (!strcmp(s,"e_ALARM"))    return e_ALARM;
    1917          21 :   if (!strcmp(s,"e_ARCH"))     return e_ARCH;
    1918          21 :   if (!strcmp(s,"e_BUG"))      return e_BUG;
    1919          21 :   if (!strcmp(s,"e_COMPONENT")) return e_COMPONENT;
    1920          21 :   if (!strcmp(s,"e_CONSTPOL")) return e_CONSTPOL;
    1921          21 :   if (!strcmp(s,"e_COPRIME"))  return e_COPRIME;
    1922          21 :   if (!strcmp(s,"e_DIM"))      return e_DIM;
    1923          21 :   if (!strcmp(s,"e_DOMAIN"))   return e_DOMAIN;
    1924          21 :   if (!strcmp(s,"e_FILE"))     return e_FILE;
    1925          21 :   if (!strcmp(s,"e_FILEDESC")) return e_FILEDESC;
    1926          21 :   if (!strcmp(s,"e_FLAG"))     return e_FLAG;
    1927          21 :   if (!strcmp(s,"e_IMPL"))     return e_IMPL;
    1928          21 :   if (!strcmp(s,"e_INV"))      return e_INV;
    1929           0 :   if (!strcmp(s,"e_IRREDPOL")) return e_IRREDPOL;
    1930           0 :   if (!strcmp(s,"e_MAXPRIME")) return e_MAXPRIME;
    1931           0 :   if (!strcmp(s,"e_MEM"))      return e_MEM;
    1932           0 :   if (!strcmp(s,"e_MISC"))     return e_MISC;
    1933           0 :   if (!strcmp(s,"e_MODULUS"))  return e_MODULUS;
    1934           0 :   if (!strcmp(s,"e_NONE"))     return e_NONE;
    1935           0 :   if (!strcmp(s,"e_NOTFUNC"))  return e_NOTFUNC;
    1936           0 :   if (!strcmp(s,"e_OP"))       return e_OP;
    1937           0 :   if (!strcmp(s,"e_OVERFLOW")) return e_OVERFLOW;
    1938           0 :   if (!strcmp(s,"e_PACKAGE"))  return e_PACKAGE;
    1939           0 :   if (!strcmp(s,"e_PREC"))     return e_PREC;
    1940           0 :   if (!strcmp(s,"e_PRIME"))    return e_PRIME;
    1941           0 :   if (!strcmp(s,"e_PRIORITY")) return e_PRIORITY;
    1942           0 :   if (!strcmp(s,"e_ROOTS0"))   return e_ROOTS0;
    1943           0 :   if (!strcmp(s,"e_SQRTN"))    return e_SQRTN;
    1944           0 :   if (!strcmp(s,"e_STACK"))    return e_STACK;
    1945           0 :   if (!strcmp(s,"e_SYNTAX"))   return e_SYNTAX;
    1946           0 :   if (!strcmp(s,"e_TYPE"))     return e_TYPE;
    1947           0 :   if (!strcmp(s,"e_TYPE2"))    return e_TYPE2;
    1948           0 :   if (!strcmp(s,"e_USER"))     return e_USER;
    1949           0 :   if (!strcmp(s,"e_VAR"))      return e_VAR;
    1950           0 :   pari_err(e_MISC,"unknown error name");
    1951             :   return -1; /* LCOV_EXCL_LINE */
    1952             : }
    1953             : 
    1954             : GEN
    1955       26927 : errname(GEN err)
    1956             : {
    1957       26927 :   if (typ(err)!=t_ERROR) pari_err_TYPE("errname",err);
    1958       26927 :   return strtoGENstr(numerr_name(err_get_num(err)));
    1959             : }
    1960             : 
    1961             : /* Try f (trapping error e), recover using r (break_loop, if NULL) */
    1962             : GEN
    1963          21 : trap0(const char *e, GEN r, GEN f)
    1964             : {
    1965          21 :   long numerr = CATCH_ALL;
    1966             :   GEN x;
    1967          21 :   if (!e || !*e) numerr = CATCH_ALL;
    1968          21 :   else numerr = name_numerr(e);
    1969          21 :   if (!f) {
    1970           0 :     pari_warn(warner,"default handlers are no longer supported --> ignored");
    1971           0 :     return gnil;
    1972             :   }
    1973          21 :   x = closure_trapgen(f, numerr);
    1974          14 :   if (x == (GEN)1L) x = r? closure_evalgen(r): gnil;
    1975          14 :   return x;
    1976             : }
    1977             : 
    1978             : /*******************************************************************/
    1979             : /*                                                                */
    1980             : /*                       CLONING & COPY                            */
    1981             : /*                  Replicate an existing GEN                      */
    1982             : /*                                                                 */
    1983             : /*******************************************************************/
    1984             : /* lontyp[tx] = 0 (non recursive type) or number of codewords for type tx */
    1985             : const  long lontyp[] = { 0,0,0,1,1,2,1,2,1,1, 2,2,0,1,1,1,1,1,1,1, 2,0,0,2,2,1 };
    1986             : 
    1987             : static GEN
    1988        1503 : list_internal_copy(GEN z, long nmax)
    1989             : {
    1990             :   long i, l;
    1991             :   GEN a;
    1992        1503 :   if (!z) return NULL;
    1993        1042 :   l = lg(z);
    1994        1042 :   a = newblock(nmax+1);
    1995       49052 :   for (i = 1; i < l; i++) gel(a,i) = gel(z,i)? gclone(gel(z,i)): gen_0;
    1996        1042 :   a[0] = z[0]; setisclone(a); return a;
    1997             : }
    1998             : 
    1999             : static void
    2000        1503 : listassign(GEN x, GEN y)
    2001             : {
    2002        1503 :   long nmax = list_nmax(x);
    2003        1503 :   GEN L = list_data(x);
    2004        1503 :   if (!nmax && L) nmax = lg(L) + 32; /* not malloc'ed yet */
    2005        1503 :   y[1] = evaltyp(list_typ(x))|evallg(nmax);
    2006        1503 :   list_data(y) = list_internal_copy(L, nmax);
    2007        1503 : }
    2008             : 
    2009             : /* transform a non-malloced list (e.g. from gtolist or gtomap) to a malloced
    2010             :  * list suitable for listput */
    2011             : GEN
    2012           0 : listinit(GEN x)
    2013             : {
    2014           0 :   GEN y = cgetg(3, t_LIST);
    2015           0 :   listassign(x, y); return y;
    2016             : }
    2017             : 
    2018             : /* copy list on the PARI stack */
    2019             : GEN
    2020         541 : listcopy(GEN x)
    2021             : {
    2022         541 :   GEN y = mklist(), L = list_data(x);
    2023         541 :   if (L) list_data(y) = gcopy(L);
    2024         541 :   y[1] = evaltyp(list_typ(x));
    2025         541 :   return y;
    2026             : }
    2027             : 
    2028             : GEN
    2029  3670141047 : gcopy(GEN x)
    2030             : {
    2031  3670141047 :   long tx = typ(x), lx, i;
    2032             :   GEN y;
    2033  3670141047 :   switch(tx)
    2034             :   { /* non recursive types */
    2035  3260017054 :     case t_INT: return signe(x)? icopy(x): gen_0;
    2036   227480058 :     case t_REAL:
    2037             :     case t_STR:
    2038   227480058 :     case t_VECSMALL: return leafcopy(x);
    2039             :     /* one more special case */
    2040         541 :     case t_LIST: return listcopy(x);
    2041             :   }
    2042   182643394 :   y = cgetg_copy(x, &lx);
    2043   182646162 :   if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2044   807244922 :   for (; i<lx; i++) gel(y,i) = gcopy(gel(x,i));
    2045   182646269 :   return y;
    2046             : }
    2047             : 
    2048             : /* as gcopy, but truncate to the first lx components if recursive type
    2049             :  * [ leaves use their own lg ]. No checks. */
    2050             : GEN
    2051         588 : gcopy_lg(GEN x, long lx)
    2052             : {
    2053         588 :   long tx = typ(x), i;
    2054             :   GEN y;
    2055         588 :   switch(tx)
    2056             :   { /* non recursive types */
    2057           0 :     case t_INT: return signe(x)? icopy(x): gen_0;
    2058           0 :     case t_REAL:
    2059             :     case t_STR:
    2060           0 :     case t_VECSMALL: return leafcopy(x);
    2061             :     /* one more special case */
    2062           0 :     case t_LIST: return listcopy(x);
    2063             :   }
    2064         588 :   y = cgetg(lx, tx);
    2065         588 :   if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2066        1694 :   for (; i<lx; i++) gel(y,i) = gcopy(gel(x,i));
    2067         588 :   return y;
    2068             : }
    2069             : 
    2070             : /* cf cgetg_copy: "allocate" (by updating first codeword only) for subsequent
    2071             :  * copy of x, as if avma = *AVMA */
    2072             : INLINE GEN
    2073   474733189 : cgetg_copy_avma(GEN x, long *plx, pari_sp *AVMA) {
    2074             :   GEN z;
    2075   474733189 :   *plx = lg(x);
    2076   474733189 :   z = ((GEN)*AVMA) - *plx;
    2077   474733189 :   z[0] = x[0] & (TYPBITS|LGBITS);
    2078   474733189 :   *AVMA = (pari_sp)z; return z;
    2079             : }
    2080             : INLINE GEN
    2081         445 : cgetlist_avma(pari_sp *AVMA)
    2082             : {
    2083         445 :   GEN y = ((GEN)*AVMA) - 3;
    2084         445 :   y[0] = _evallg(3) | evaltyp(t_LIST);
    2085         445 :   *AVMA = (pari_sp)y; return y;
    2086             : }
    2087             : 
    2088             : /* copy x as if avma = *AVMA, update *AVMA */
    2089             : GEN
    2090  2843927756 : gcopy_avma(GEN x, pari_sp *AVMA)
    2091             : {
    2092  2843927756 :   long i, lx, tx = typ(x);
    2093             :   GEN y;
    2094             : 
    2095  2843927756 :   switch(typ(x))
    2096             :   { /* non recursive types */
    2097  2698464881 :     case t_INT:
    2098  2698464881 :       if (lgefint(x) == 2) return gen_0;
    2099  2213227159 :       *AVMA = (pari_sp)icopy_avma(x, *AVMA);
    2100  2213227220 :       return (GEN)*AVMA;
    2101    33150423 :     case t_REAL: case t_STR: case t_VECSMALL:
    2102    33150423 :       *AVMA = (pari_sp)leafcopy_avma(x, *AVMA);
    2103    33150423 :       return (GEN)*AVMA;
    2104             : 
    2105             :     /* one more special case */
    2106         441 :     case t_LIST:
    2107         441 :       y = cgetlist_avma(AVMA);
    2108         441 :       listassign(x, y); return y;
    2109             : 
    2110             :   }
    2111   112312011 :   y = cgetg_copy_avma(x, &lx, AVMA);
    2112   112312032 :   if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2113   519314234 :   for (; i<lx; i++) gel(y,i) = gcopy_avma(gel(x,i), AVMA);
    2114   112312035 :   return y;
    2115             : }
    2116             : 
    2117             : /* [copy_bin/bin_copy:] same as gcopy_avma but use NULL to code an exact 0, and
    2118             :  * make shallow copies of finalized t_LISTs */
    2119             : static GEN
    2120  2437029735 : gcopy_av0(GEN x, pari_sp *AVMA)
    2121             : {
    2122  2437029735 :   long i, lx, tx = typ(x);
    2123             :   GEN y;
    2124             : 
    2125  2437029735 :   switch(tx)
    2126             :   { /* non recursive types */
    2127  1858053806 :     case t_INT:
    2128  1858053806 :       if (!signe(x)) return NULL; /* special marker */
    2129   889497100 :       *AVMA = (pari_sp)icopy_avma(x, *AVMA);
    2130   889502761 :       return (GEN)*AVMA;
    2131          49 :     case t_LIST:
    2132          49 :       if (list_data(x) && !list_nmax(x)) break; /* not finalized, need copy */
    2133             :       /* else finalized: shallow copy */
    2134             :     case t_REAL: case t_STR: case t_VECSMALL:
    2135   216561212 :       *AVMA = (pari_sp)leafcopy_avma(x, *AVMA);
    2136   216559604 :       return (GEN)*AVMA;
    2137             :   }
    2138   362414717 :   y = cgetg_copy_avma(x, &lx, AVMA);
    2139   362421145 :   if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2140  2503053853 :   for (; i<lx; i++) gel(y,i) = gcopy_av0(gel(x,i), AVMA);
    2141   362422458 :   return y;
    2142             : }
    2143             : 
    2144             : INLINE GEN
    2145          12 : icopy_avma_canon(GEN x, pari_sp AVMA)
    2146             : {
    2147          12 :   long i, lx = lgefint(x);
    2148          12 :   GEN y = ((GEN)AVMA) - lx;
    2149          12 :   y[0] = evaltyp(t_INT)|evallg(lx); /* kills isclone */
    2150          12 :   y[1] = x[1]; x = int_MSW(x);
    2151          24 :   for (i=2; i<lx; i++, x = int_precW(x)) y[i] = *x;
    2152          12 :   return y;
    2153             : }
    2154             : 
    2155             : /* [copy_bin_canon:] same as gcopy_av0, but copy integers in
    2156             :  * canonical (native kernel) form and make a full copy of t_LISTs */
    2157             : static GEN
    2158          32 : gcopy_av0_canon(GEN x, pari_sp *AVMA)
    2159             : {
    2160          32 :   long i, lx, tx = typ(x);
    2161             :   GEN y;
    2162             : 
    2163          32 :   switch(tx)
    2164             :   { /* non recursive types */
    2165          20 :     case t_INT:
    2166          20 :       if (!signe(x)) return NULL; /* special marker */
    2167          12 :       *AVMA = (pari_sp)icopy_avma_canon(x, *AVMA);
    2168          12 :       return (GEN)*AVMA;
    2169           0 :     case t_REAL: case t_STR: case t_VECSMALL:
    2170           0 :       *AVMA = (pari_sp)leafcopy_avma(x, *AVMA);
    2171           0 :       return (GEN)*AVMA;
    2172             : 
    2173             :     /* one more special case */
    2174           4 :     case t_LIST:
    2175             :     {
    2176           4 :       long t = list_typ(x);
    2177           4 :       GEN y = cgetlist_avma(AVMA), z = list_data(x);
    2178           4 :       if (z) {
    2179           0 :         list_data(y) = gcopy_av0_canon(z, AVMA);
    2180           0 :         y[1] = evaltyp(t)|evallg(lg(z)-1);
    2181             :       } else {
    2182           4 :         list_data(y) = NULL;
    2183           4 :         y[1] = evaltyp(t);
    2184             :       }
    2185           4 :       return y;
    2186             :     }
    2187             :   }
    2188           8 :   y = cgetg_copy_avma(x, &lx, AVMA);
    2189           8 :   if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2190          24 :   for (; i<lx; i++) gel(y,i) = gcopy_av0_canon(gel(x,i), AVMA);
    2191           8 :   return y;
    2192             : }
    2193             : 
    2194             : /* [copy_bin/bin_copy:] size (number of words) required for
    2195             :  * gcopy_av0_canon(x) */
    2196             : static long
    2197          32 : taille0_canon(GEN x)
    2198             : {
    2199          32 :   long i,n,lx, tx = typ(x);
    2200          32 :   switch(tx)
    2201             :   { /* non recursive types */
    2202          20 :     case t_INT: return signe(x)? lgefint(x): 0;
    2203           0 :     case t_REAL:
    2204             :     case t_STR:
    2205           0 :     case t_VECSMALL: return lg(x);
    2206             : 
    2207             :     /* one more special case */
    2208           4 :     case t_LIST:
    2209             :     {
    2210           4 :       GEN L = list_data(x);
    2211           4 :       return L? 3 + taille0_canon(L): 3;
    2212             :     }
    2213             :   }
    2214           8 :   n = lx = lg(x);
    2215          24 :   for (i=lontyp[tx]; i<lx; i++) n += taille0_canon(gel(x,i));
    2216           8 :   return n;
    2217             : }
    2218             : 
    2219             : /* [copy_bin/bin_copy:] size (number of words) required for gcopy_av0(x) */
    2220             : static long
    2221  2437061671 : taille0(GEN x)
    2222             : {
    2223  2437061671 :   long i,n,lx, tx = typ(x);
    2224  2437061671 :   switch(tx)
    2225             :   { /* non recursive types */
    2226  1858067671 :     case t_INT:
    2227  1858067671 :       lx = lgefint(x);
    2228  1858067671 :       return lx == 2? 0: lx;
    2229          49 :     case t_LIST:
    2230             :     {
    2231          49 :       GEN L = list_data(x);
    2232          49 :       if (L && !list_nmax(x)) break; /* not finalized, deep copy */
    2233             :     }
    2234             :     /* else finalized: shallow */
    2235             :     case t_REAL:
    2236             :     case t_STR:
    2237             :     case t_VECSMALL:
    2238   216580580 :       return lg(x);
    2239             :   }
    2240   362413420 :   n = lx = lg(x);
    2241  2503080524 :   for (i=lontyp[tx]; i<lx; i++) n += taille0(gel(x,i));
    2242   362416595 :   return n;
    2243             : }
    2244             : 
    2245             : static long
    2246  2917087202 : gsizeclone_i(GEN x)
    2247             : {
    2248  2917087202 :   long i,n,lx, tx = typ(x);
    2249  2917087202 :   switch(tx)
    2250             :   { /* non recursive types */
    2251  2699039434 :     case t_INT: lx = lgefint(x); return lx == 2? 0: lx;;
    2252    40837038 :     case t_REAL:
    2253             :     case t_STR:
    2254    40837038 :     case t_VECSMALL: return lg(x);
    2255             : 
    2256        1503 :     case t_LIST: return 3;
    2257   177209227 :     default:
    2258   177209227 :       n = lx = lg(x);
    2259  3023351126 :       for (i=lontyp[tx]; i<lx; i++) n += gsizeclone_i(gel(x,i));
    2260   177209303 :       return n;
    2261             :   }
    2262             : }
    2263             : 
    2264             : /* #words needed to clone x; t_LIST is a special case since list_data() is
    2265             :  * malloc'ed later, in list_internal_copy() */
    2266             : static long
    2267   211307629 : gsizeclone(GEN x) { return (typ(x) == t_INT)? lgefint(x): gsizeclone_i(x); }
    2268             : 
    2269             : long
    2270         147 : gsizeword(GEN x)
    2271             : {
    2272             :   GEN L;
    2273         147 :   if (typ(x) != t_LIST) return gsizeclone(x);
    2274             :   /* For t_LIST, return the actual list size, gsizeclone() is always 3 */
    2275           0 :   L = list_data(x);
    2276           0 :   return L? 3 + gsizeclone(L): 3;
    2277             : }
    2278             : long
    2279         147 : gsizebyte(GEN x) { return gsizeword(x) * sizeof(long); }
    2280             : 
    2281             : /* return a clone of x structured as a gcopy */
    2282             : GENbin*
    2283   296397417 : copy_bin(GEN x)
    2284             : {
    2285   296397417 :   long t = taille0(x);
    2286   296398473 :   GENbin *p = (GENbin*)pari_malloc(sizeof(GENbin) + t*sizeof(long));
    2287   296399322 :   pari_sp AVMA = (pari_sp)(GENbinbase(p) + t);
    2288   296399066 :   p->rebase = &shiftaddress;
    2289   296399066 :   p->len = t;
    2290   296399066 :   p->x   = gcopy_av0(x, &AVMA);
    2291   296398716 :   p->base= (GEN)AVMA; return p;
    2292             : }
    2293             : 
    2294             : /* same, writing t_INT in canonical native form */
    2295             : GENbin*
    2296          16 : copy_bin_canon(GEN x)
    2297             : {
    2298          16 :   long t = taille0_canon(x);
    2299          16 :   GENbin *p = (GENbin*)pari_malloc(sizeof(GENbin) + t*sizeof(long));
    2300          16 :   pari_sp AVMA = (pari_sp)(GENbinbase(p) + t);
    2301          16 :   p->rebase = &shiftaddress_canon;
    2302          16 :   p->len = t;
    2303          16 :   p->x   = gcopy_av0_canon(x, &AVMA);
    2304          16 :   p->base= (GEN)AVMA; return p;
    2305             : }
    2306             : 
    2307             : GEN
    2308   211307487 : gclone(GEN x)
    2309             : {
    2310   211307487 :   long i,lx,tx = typ(x), t = gsizeclone(x);
    2311   211307609 :   GEN y = newblock(t);
    2312   211307477 :   switch(tx)
    2313             :   { /* non recursive types */
    2314   140362089 :     case t_INT:
    2315   140362089 :       lx = lgefint(x);
    2316   140362089 :       y[0] = evaltyp(t_INT)|evallg(lx);
    2317   759298055 :       for (i=1; i<lx; i++) y[i] = x[i];
    2318   140362076 :       break;
    2319     6584987 :     case t_REAL:
    2320             :     case t_STR:
    2321             :     case t_VECSMALL:
    2322     6584987 :       lx = lg(x);
    2323    82566853 :       for (i=0; i<lx; i++) y[i] = x[i];
    2324     6584987 :       break;
    2325             : 
    2326             :     /* one more special case */
    2327        1062 :     case t_LIST:
    2328        1062 :       y[0] = evaltyp(t_LIST)|_evallg(3);
    2329        1062 :       listassign(x, y);
    2330        1062 :       break;
    2331    64359339 :     default: {
    2332    64359339 :       pari_sp AVMA = (pari_sp)(y + t);
    2333    64359339 :       lx = lg(x);
    2334    64359339 :       y[0] = x[0];
    2335    64359339 :       if (lontyp[tx] == 1) i = 1; else { y[1] = x[1]; i = 2; }
    2336  2501284939 :       for (; i<lx; i++) gel(y,i) = gcopy_avma(gel(x,i), &AVMA);
    2337             :     }
    2338             :   }
    2339   211307486 :   setisclone(y); return y;
    2340             : }
    2341             : 
    2342             : void
    2343  1468508705 : shiftaddress(GEN x, long dec)
    2344             : {
    2345  1468508705 :   long i, lx, tx = typ(x);
    2346  1468508705 :   if (is_recursive_t(tx))
    2347             :   {
    2348   362421720 :     if (tx == t_LIST)
    2349             :     {
    2350          49 :       if (!list_data(x) || list_nmax(x)) return; /* empty or finalized */
    2351             :       /* not finalized, update pointers  */
    2352             :     }
    2353   362421678 :     lx = lg(x);
    2354  2503094504 :     for (i=lontyp[tx]; i<lx; i++) {
    2355  2140673739 :       if (!x[i]) gel(x,i) = gen_0;
    2356             :       else
    2357             :       {
    2358  1196014487 :         x[i] += dec;
    2359  1196014487 :         shiftaddress(gel(x,i), dec);
    2360             :       }
    2361             :     }
    2362             :   }
    2363             : }
    2364             : 
    2365             : void
    2366          24 : shiftaddress_canon(GEN x, long dec)
    2367             : {
    2368          24 :   long i, lx, tx = typ(x);
    2369          24 :   switch(tx)
    2370             :   { /* non recursive types */
    2371          12 :     case t_INT: {
    2372             :       GEN y;
    2373          12 :       lx = lgefint(x); if (lx <= 3) return;
    2374           0 :       y = x + 2;
    2375           0 :       x = int_MSW(x);  if (x == y) return;
    2376           0 :       while (x > y) { lswap(*x, *y); x = int_precW(x); y++; }
    2377           0 :       break;
    2378             :     }
    2379           0 :     case t_REAL:
    2380             :     case t_STR:
    2381             :     case t_VECSMALL:
    2382           0 :       break;
    2383             : 
    2384             :     /* one more special case */
    2385           4 :     case t_LIST:
    2386           4 :       if (!list_data(x)) break;
    2387             :     default: /* Fall through */
    2388           8 :       lx = lg(x);
    2389          24 :       for (i=lontyp[tx]; i<lx; i++) {
    2390          16 :         if (!x[i]) gel(x,i) = gen_0;
    2391             :         else
    2392             :         {
    2393           8 :           x[i] += dec;
    2394           8 :           shiftaddress_canon(gel(x,i), dec);
    2395             :         }
    2396             :       }
    2397             :   }
    2398             : }
    2399             : 
    2400             : /********************************************************************/
    2401             : /**                                                                **/
    2402             : /**                INSERT DYNAMIC OBJECT IN STRUCTURE              **/
    2403             : /**                                                                **/
    2404             : /********************************************************************/
    2405             : GEN
    2406          35 : obj_reinit(GEN S)
    2407             : {
    2408          35 :   GEN s, T = leafcopy(S);
    2409          35 :   long a = lg(T)-1;
    2410          35 :   s = gel(T,a);
    2411          35 :   gel(T,a) = zerovec(lg(s)-1);
    2412          35 :   return T;
    2413             : }
    2414             : 
    2415             : GEN
    2416     1217084 : obj_init(long d, long n)
    2417             : {
    2418     1217084 :   GEN S = cgetg(d+2, t_VEC);
    2419     1217084 :   gel(S, d+1) = zerovec(n);
    2420     1217084 :   return S;
    2421             : }
    2422             : /* insert O in S [last position] at position K, return it */
    2423             : GEN
    2424     1127882 : obj_insert(GEN S, long K, GEN O)
    2425     1127882 : { return obj_insert_shallow(S, K, gclone(O)); }
    2426             : /* as obj_insert. WITHOUT cloning (for libpari, when creating a *new* obj
    2427             :  * from an existing one) */
    2428             : GEN
    2429     1134406 : obj_insert_shallow(GEN S, long K, GEN O)
    2430             : {
    2431     1134406 :   GEN o, v = gel(S, lg(S)-1);
    2432     1134406 :   if (typ(v) != t_VEC) pari_err_TYPE("obj_insert", S);
    2433     1134406 :   o = gel(v,K);
    2434     1134406 :   gel(v,K) = O; /*SIGINT: before unclone(o)*/
    2435     1134406 :   if (isclone(o)) gunclone(o); return gel(v,K);
    2436             : }
    2437             : 
    2438             : /* Does S [last position] contain data at position K ? Return it, or NULL */
    2439             : GEN
    2440     2718363 : obj_check(GEN S, long K)
    2441             : {
    2442     2718363 :   GEN O, v = gel(S, lg(S)-1);
    2443     2718363 :   if (typ(v) != t_VEC || K >= lg(v)) pari_err_TYPE("obj_check", S);
    2444     2718363 :   O = gel(v,K); return isintzero(O)? NULL: O;
    2445             : }
    2446             : 
    2447             : GEN
    2448      863899 : obj_checkbuild(GEN S, long tag, GEN (*build)(GEN))
    2449             : {
    2450      863899 :   GEN O = obj_check(S, tag);
    2451      863899 :   if (!O)
    2452      657105 :   { pari_sp av = avma; O = obj_insert(S, tag, build(S)); set_avma(av); }
    2453      863892 :   return O;
    2454             : }
    2455             : 
    2456             : GEN
    2457      107999 : obj_checkbuild_prec(GEN S, long tag, GEN (*build)(GEN,long),
    2458             :   long (*pr)(GEN), long prec)
    2459             : {
    2460      107999 :   pari_sp av = avma;
    2461      107999 :   GEN w = obj_check(S, tag);
    2462      107999 :   if (!w || pr(w) < prec) w = obj_insert(S, tag, build(S, prec));
    2463      107999 :   set_avma(av); return gcopy(w);
    2464             : }
    2465             : GEN
    2466       13667 : obj_checkbuild_realprec(GEN S, long tag, GEN (*build)(GEN,long), long prec)
    2467       13667 : { return obj_checkbuild_prec(S,tag,build,gprecision,prec); }
    2468             : GEN
    2469         497 : obj_checkbuild_padicprec(GEN S, long tag, GEN (*build)(GEN,long), long prec)
    2470         497 : { return obj_checkbuild_prec(S,tag,build,padicprec_relative,prec); }
    2471             : 
    2472             : /* Reset S [last position], freeing all clones */
    2473             : void
    2474       13265 : obj_free(GEN S)
    2475             : {
    2476       13265 :   GEN v = gel(S, lg(S)-1);
    2477             :   long i;
    2478       13265 :   if (typ(v) != t_VEC) pari_err_TYPE("obj_free", S);
    2479       74522 :   for (i = 1; i < lg(v); i++)
    2480             :   {
    2481       61257 :     GEN o = gel(v,i);
    2482       61257 :     gel(v,i) = gen_0;
    2483       61257 :     gunclone_deep(o);
    2484             :   }
    2485       13265 : }
    2486             : 
    2487             : /*******************************************************************/
    2488             : /*                                                                 */
    2489             : /*                         STACK MANAGEMENT                        */
    2490             : /*                                                                 */
    2491             : /*******************************************************************/
    2492             : INLINE void
    2493  2646258651 : dec_gerepile(pari_sp *x, pari_sp av0, pari_sp av, pari_sp tetpil, size_t dec)
    2494             : {
    2495  2646258651 :   if (*x < av && *x >= av0)
    2496             :   { /* update address if in stack */
    2497  2229654500 :     if (*x < tetpil) *x += dec;
    2498           0 :     else pari_err_BUG("gerepile, significant pointers lost");
    2499             :   }
    2500  2646259447 : }
    2501             : 
    2502             : void
    2503     1165932 : gerepileallsp(pari_sp av, pari_sp tetpil, int n, ...)
    2504             : {
    2505     1165932 :   const pari_sp av0 = avma;
    2506     1165932 :   const size_t dec = av-tetpil;
    2507             :   int i;
    2508     1165932 :   va_list a; va_start(a, n);
    2509     1165932 :   (void)gerepile(av,tetpil,NULL);
    2510     5342504 :   for (i=0; i<n; i++) dec_gerepile((pari_sp*)va_arg(a,GEN*), av0,av,tetpil,dec);
    2511     1165932 :   va_end(a);
    2512     1165932 : }
    2513             : 
    2514             : /* Takes an array of pointers to GENs, of length n.
    2515             :  * Cleans up the stack between av and tetpil, updating those GENs. */
    2516             : void
    2517     6460730 : gerepilemanysp(pari_sp av, pari_sp tetpil, GEN* gptr[], int n)
    2518             : {
    2519     6460730 :   const pari_sp av0 = avma;
    2520     6460730 :   const size_t dec = av-tetpil;
    2521             :   int i;
    2522     6460730 :   (void)gerepile(av,tetpil,NULL);
    2523    19385564 :   for (i=0; i<n; i++) dec_gerepile((pari_sp*)gptr[i], av0, av, tetpil, dec);
    2524     6460730 : }
    2525             : 
    2526             : /* Takes an array of GENs (cast to longs), of length n.
    2527             :  * Cleans up the stack between av and tetpil, updating those GENs. */
    2528             : void
    2529   154095262 : gerepilecoeffssp(pari_sp av, pari_sp tetpil, long *g, int n)
    2530             : {
    2531   154095262 :   const pari_sp av0 = avma;
    2532   154095262 :   const size_t dec = av-tetpil;
    2533             :   int i;
    2534   154095262 :   (void)gerepile(av,tetpil,NULL);
    2535   462285786 :   for (i=0; i<n; i++,g++) dec_gerepile((pari_sp*)g, av0, av, tetpil, dec);
    2536   154095262 : }
    2537             : 
    2538             : static int
    2539           0 : dochk_gerepileupto(GEN av, GEN x)
    2540             : {
    2541             :   long i,lx,tx;
    2542           0 :   if (!isonstack(x)) return 1;
    2543           0 :   if (x > av)
    2544             :   {
    2545           0 :     pari_warn(warner,"bad object %Ps",x);
    2546           0 :     return 0;
    2547             :   }
    2548           0 :   tx = typ(x);
    2549           0 :   if (! is_recursive_t(tx)) return 1;
    2550             : 
    2551           0 :   lx = lg(x);
    2552           0 :   for (i=lontyp[tx]; i<lx; i++)
    2553           0 :     if (!dochk_gerepileupto(av, gel(x,i)))
    2554             :     {
    2555           0 :       pari_warn(warner,"bad component %ld in object %Ps",i,x);
    2556           0 :       return 0;
    2557             :     }
    2558           0 :   return 1;
    2559             : }
    2560             : /* check that x and all its components are out of stack, or have been
    2561             :  * created after av */
    2562             : int
    2563           0 : chk_gerepileupto(GEN x) { return dochk_gerepileupto(x, x); }
    2564             : 
    2565             : /* print stack between avma & av */
    2566             : void
    2567           0 : dbg_gerepile(pari_sp av)
    2568             : {
    2569           0 :   GEN x = (GEN)avma;
    2570           0 :   while (x < (GEN)av)
    2571             :   {
    2572           0 :     const long tx = typ(x), lx = lg(x);
    2573             :     GEN *a;
    2574             : 
    2575           0 :     pari_printf(" [%ld] %Ps:", x - (GEN)avma, x);
    2576           0 :     if (! is_recursive_t(tx)) { pari_putc('\n'); x += lx; continue; }
    2577           0 :     a = (GEN*)x + lontyp[tx]; x += lx;
    2578           0 :     for (  ; a < (GEN*)x; a++)
    2579             :     {
    2580           0 :       if (*a == gen_0)
    2581           0 :         pari_puts("  gen_0");
    2582           0 :       else if (*a == gen_1)
    2583           0 :         pari_puts("  gen_1");
    2584           0 :       else if (*a == gen_m1)
    2585           0 :         pari_puts("  gen_m1");
    2586           0 :       else if (*a == gen_2)
    2587           0 :         pari_puts("  gen_2");
    2588           0 :       else if (*a == gen_m2)
    2589           0 :         pari_puts("  gen_m2");
    2590           0 :       else if (*a == ghalf)
    2591           0 :         pari_puts("  ghalf");
    2592           0 :       else if (isclone(*a))
    2593           0 :         pari_printf("  %Ps (clone)", *a);
    2594             :       else
    2595           0 :         pari_printf("  %Ps [%ld]", *a, *a - (GEN)avma);
    2596           0 :       if (a+1 < (GEN*)x) pari_putc(',');
    2597             :     }
    2598           0 :     pari_printf("\n");
    2599             :   }
    2600           0 : }
    2601             : void
    2602           0 : dbg_gerepileupto(GEN q)
    2603             : {
    2604           0 :   err_printf("%Ps:\n", q);
    2605           0 :   dbg_gerepile((pari_sp) (q+lg(q)));
    2606           0 : }
    2607             : 
    2608             : GEN
    2609   580379922 : gerepile(pari_sp av, pari_sp tetpil, GEN q)
    2610             : {
    2611   580379922 :   const size_t dec = av - tetpil;
    2612   580379922 :   const pari_sp av0 = avma;
    2613             :   GEN x, a;
    2614             : 
    2615   580379922 :   if (dec == 0) return q;
    2616   482905267 :   if ((long)dec < 0) pari_err(e_MISC,"lbot>ltop in gerepile");
    2617             : 
    2618             :   /* dec_gerepile(&q, av0, av, tetpil, dec), saving 1 comparison */
    2619   482908339 :   if (q >= (GEN)av0 && q < (GEN)tetpil)
    2620   325087734 :     q = (GEN) (((pari_sp)q) + dec);
    2621             : 
    2622 14737664555 :   for (x = (GEN)av, a = (GEN)tetpil; a > (GEN)av0; ) *--x = *--a;
    2623   482908339 :   set_avma((pari_sp)x);
    2624  3104742730 :   while (x < (GEN)av)
    2625             :   {
    2626  2621825840 :     const long tx = typ(x), lx = lg(x);
    2627             : 
    2628  2621825840 :     if (! is_recursive_t(tx)) { x += lx; continue; }
    2629   570509080 :     a = x + lontyp[tx]; x += lx;
    2630  2891473968 :     for (  ; a < x; a++) dec_gerepile((pari_sp*)a, av0, av, tetpil, dec);
    2631             :   }
    2632   482916890 :   return q;
    2633             : }
    2634             : 
    2635             : void
    2636           0 : fill_stack(void)
    2637             : {
    2638           0 :   GEN x = ((GEN)pari_mainstack->bot);
    2639           0 :   while (x < (GEN)avma) *x++ = 0xfefefefeUL;
    2640           0 : }
    2641             : 
    2642             : void
    2643           0 : debug_stack(void)
    2644             : {
    2645           0 :   pari_sp top = pari_mainstack->top, bot = pari_mainstack->bot;
    2646             :   GEN z;
    2647           0 :   err_printf("bot=0x%lx\ttop=0x%lx\tavma=0x%lx\n", bot, top, avma);
    2648           0 :   for (z = ((GEN)top)-1; z >= (GEN)avma; z--)
    2649           0 :     err_printf("%p:\t0x%lx\t%lu\n",z,*z,*z);
    2650           0 : }
    2651             : 
    2652             : void
    2653           0 : setdebugvar(long n) { DEBUGVAR=n; }
    2654             : 
    2655             : long
    2656           0 : getdebugvar(void) { return DEBUGVAR; }
    2657             : 
    2658             : long
    2659           7 : getstack(void) { return pari_mainstack->top-avma; }
    2660             : 
    2661             : /*******************************************************************/
    2662             : /*                                                                 */
    2663             : /*                               timer_delay                             */
    2664             : /*                                                                 */
    2665             : /*******************************************************************/
    2666             : 
    2667             : #if defined(USE_CLOCK_GETTIME)
    2668             : #if defined(_POSIX_THREAD_CPUTIME)
    2669             : static THREAD clockid_t time_type = CLOCK_THREAD_CPUTIME_ID;
    2670             : #else
    2671             : static const THREAD clockid_t time_type = CLOCK_PROCESS_CPUTIME_ID;
    2672             : #endif
    2673             : static void
    2674             : pari_init_timer(void)
    2675             : {
    2676             : #if defined(_POSIX_THREAD_CPUTIME)
    2677             :   time_type = CLOCK_PROCESS_CPUTIME_ID;
    2678             : #endif
    2679             : }
    2680             : 
    2681             : void
    2682             : timer_start(pari_timer *T)
    2683             : {
    2684             :   struct timespec t;
    2685             :   clock_gettime(time_type,&t);
    2686             :   T->us = t.tv_nsec / 1000;
    2687             :   T->s  = t.tv_sec;
    2688             : }
    2689             : #elif defined(USE_GETRUSAGE)
    2690             : #ifdef RUSAGE_THREAD
    2691             : static THREAD int rusage_type = RUSAGE_THREAD;
    2692             : #else
    2693             : static const THREAD int rusage_type = RUSAGE_SELF;
    2694             : #endif /*RUSAGE_THREAD*/
    2695             : static void
    2696        1690 : pari_init_timer(void)
    2697             : {
    2698             : #ifdef RUSAGE_THREAD
    2699        1690 :   rusage_type = RUSAGE_SELF;
    2700             : #endif
    2701        1690 : }
    2702             : 
    2703             : void
    2704      256181 : timer_start(pari_timer *T)
    2705             : {
    2706             :   struct rusage r;
    2707      256181 :   getrusage(rusage_type,&r);
    2708      256181 :   T->us = r.ru_utime.tv_usec;
    2709      256181 :   T->s  = r.ru_utime.tv_sec;
    2710      256181 : }
    2711             : #elif defined(USE_FTIME)
    2712             : 
    2713             : static void
    2714             : pari_init_timer(void) { }
    2715             : 
    2716             : void
    2717             : timer_start(pari_timer *T)
    2718             : {
    2719             :   struct timeb t;
    2720             :   ftime(&t);
    2721             :   T->us = ((long)t.millitm) * 1000;
    2722             :   T->s  = t.time;
    2723             : }
    2724             : 
    2725             : #else
    2726             : 
    2727             : static void
    2728             : _get_time(pari_timer *T, long Ticks, long TickPerSecond)
    2729             : {
    2730             :   T->us = (long) ((Ticks % TickPerSecond) * (1000000. / TickPerSecond));
    2731             :   T->s  = Ticks / TickPerSecond;
    2732             : }
    2733             : 
    2734             : # ifdef USE_TIMES
    2735             : static void
    2736             : pari_init_timer(void) { }
    2737             : 
    2738             : void
    2739             : timer_start(pari_timer *T)
    2740             : {
    2741             : # ifdef _SC_CLK_TCK
    2742             :   long tck = sysconf(_SC_CLK_TCK);
    2743             : # else
    2744             :   long tck = CLK_TCK;
    2745             : # endif
    2746             :   struct tms t; times(&t);
    2747             :   _get_time(T, t.tms_utime, tck);
    2748             : }
    2749             : # elif defined(_WIN32)
    2750             : static void
    2751             : pari_init_timer(void) { }
    2752             : 
    2753             : void
    2754             : timer_start(pari_timer *T)
    2755             : { _get_time(T, win32_timer(), 1000); }
    2756             : # else
    2757             : #  include <time.h>
    2758             : #  ifndef CLOCKS_PER_SEC
    2759             : #   define CLOCKS_PER_SEC 1000000 /* may be false on YOUR system */
    2760             : #  endif
    2761             : static void
    2762             : pari_init_timer(void) { }
    2763             : 
    2764             : void
    2765             : timer_start(pari_timer *T)
    2766             : { _get_time(T, clock(), CLOCKS_PER_SEC); }
    2767             : # endif
    2768             : #endif
    2769             : 
    2770             : /* round microseconds to milliseconds */
    2771             : static long
    2772      164690 : rndus(long x) { return (x + 500) / 1000; }
    2773             : static long
    2774      164682 : timer_aux(pari_timer *T, pari_timer *U, void (*settime)(pari_timer *))
    2775             : {
    2776      164682 :   long s = T->s, us = T->us;
    2777      164682 :   settime(U); return 1000 * (U->s - s) + rndus(U->us - us);
    2778             : }
    2779             : 
    2780             : /* return delay, set timer checkpoint */
    2781             : long
    2782       83136 : timer_delay(pari_timer *T) { return timer_aux(T, T, &timer_start); }
    2783             : /* return delay, don't set checkpoint */
    2784             : long
    2785        1708 : timer_get(pari_timer *T) {pari_timer t; return timer_aux(T, &t, &timer_start);}
    2786             : 
    2787             : static void
    2788           0 : timer_vprintf(pari_timer *T, const char *format, va_list args)
    2789             : {
    2790           0 :   out_puts(pariErr, "Time ");
    2791           0 :   out_vprintf(pariErr, format,args);
    2792           0 :   out_printf(pariErr, ": %ld\n", timer_delay(T));
    2793           0 :   pariErr->flush();
    2794           0 : }
    2795             : void
    2796           0 : timer_printf(pari_timer *T, const char *format, ...)
    2797             : {
    2798           0 :   va_list args; va_start(args, format);
    2799           0 :   timer_vprintf(T, format, args);
    2800           0 :   va_end(args);
    2801           0 : }
    2802             : 
    2803             : long
    2804           0 : timer(void)  { static THREAD pari_timer T; return timer_delay(&T);}
    2805             : long
    2806        3298 : gettime(void)  { static THREAD pari_timer T; return timer_delay(&T);}
    2807             : 
    2808             : static THREAD pari_timer timer2_T, abstimer_T;
    2809             : long
    2810           0 : timer2(void) {  return timer_delay(&timer2_T);}
    2811             : void
    2812           0 : msgtimer(const char *format, ...)
    2813             : {
    2814           0 :   va_list args; va_start(args, format);
    2815           0 :   timer_vprintf(&timer2_T, format, args);
    2816           0 :   va_end(args);
    2817           0 : }
    2818             : long
    2819        1702 : getabstime(void)  { return timer_get(&abstimer_T);}
    2820             : 
    2821             : void
    2822      188197 : walltimer_start(pari_timer *ti)
    2823             : {
    2824             : #if defined(USE_CLOCK_GETTIME)
    2825             :   struct timespec t;
    2826             :   if (!clock_gettime(CLOCK_REALTIME,&t))
    2827             :   { ti->s = t.tv_sec; ti->us = rndus(t.tv_nsec); return; }
    2828             : #elif defined(USE_GETTIMEOFDAY)
    2829             :   struct timeval tv;
    2830      188197 :   if (!gettimeofday(&tv, NULL))
    2831      188197 :   {  ti->s = tv.tv_sec; ti->us = tv.tv_usec; return; }
    2832             : #elif defined(USE_FTIMEFORWALLTIME)
    2833             :   struct timeb tp;
    2834             :   if (!ftime(&tp))
    2835             :   { ti->s = tp.time; ti->us = tp.millitm*1000; return; }
    2836             : #endif
    2837           0 :   timer_start(ti);
    2838             : }
    2839             : /* return delay, set timer checkpoint */
    2840             : long
    2841       79838 : walltimer_delay(pari_timer *T) { return timer_aux(T, T, &walltimer_start); }
    2842             : /* return delay, don't set checkpoint */
    2843             : long
    2844           0 : walltimer_get(pari_timer *T)
    2845             : {
    2846             :   pari_timer t;
    2847           0 :   return timer_aux(T, &t, &walltimer_start);
    2848             : }
    2849             : 
    2850             : static GEN
    2851           8 : timetoi(ulong s, ulong m)
    2852             : {
    2853           8 :   pari_sp av = avma;
    2854           8 :   return gerepileuptoint(av, addiu(muluu(s, 1000), m));
    2855             : }
    2856             : GEN
    2857           8 : getwalltime(void)
    2858             : {
    2859             :   pari_timer ti;
    2860           8 :   walltimer_start(&ti);
    2861           8 :   return timetoi(ti.s, rndus(ti.us));
    2862             : }
    2863             : 
    2864             : /*******************************************************************/
    2865             : /*                                                                 */
    2866             : /*                   FUNCTIONS KNOWN TO THE ANALYZER               */
    2867             : /*                                                                 */
    2868             : /*******************************************************************/
    2869             : GEN
    2870           7 : pari_version(void)
    2871             : {
    2872           7 :   const ulong mask = (1UL<<PARI_VERSION_SHIFT) - 1;
    2873           7 :   ulong major, minor, patch, n = paricfg_version_code;
    2874           7 :   patch = n & mask; n >>= PARI_VERSION_SHIFT;
    2875           7 :   minor = n & mask; n >>= PARI_VERSION_SHIFT;
    2876           7 :   major = n;
    2877           7 :   if (*paricfg_vcsversion) {
    2878           7 :     const char *ver = paricfg_vcsversion;
    2879           7 :     const char *s = strchr(ver, '-');
    2880             :     char t[8];
    2881           7 :     const long len = s-ver;
    2882             :     GEN v;
    2883           7 :     if (!s || len > 6) pari_err_BUG("pari_version()"); /* paranoia */
    2884           7 :     memcpy(t, ver, len); t[len] = 0;
    2885           7 :     v = cgetg(6, t_VEC);
    2886           7 :     gel(v,1) = utoi(major);
    2887           7 :     gel(v,2) = utoi(minor);
    2888           7 :     gel(v,3) = utoi(patch);
    2889           7 :     gel(v,4) = stoi( atoi(t) );
    2890           7 :     gel(v,5) = strtoGENstr(s+1);
    2891           7 :     return v;
    2892             :   } else {
    2893           0 :     GEN v = cgetg(4, t_VEC);
    2894           0 :     gel(v,1) = utoi(major);
    2895           0 :     gel(v,2) = utoi(minor);
    2896           0 :     gel(v,3) = utoi(patch);
    2897           0 :     return v;
    2898             :   }
    2899             : }
    2900             : 
    2901             : /* List of GP functions: generated from the description system.
    2902             :  * Format (struct entree) :
    2903             :  *   char *name   : name (under GP).
    2904             :  *   ulong valence: (EpNEW, EpALIAS,EpVAR, EpINSTALL)|EpSTATIC
    2905             :  *   void *value  : For PREDEFINED FUNCTIONS: C function to call.
    2906             :  *                  For USER FUNCTIONS: pointer to defining data (block) =
    2907             :  *                   entree*: NULL, list of entree (arguments), NULL
    2908             :  *                   char*  : function text
    2909             :  *   long menu    : which help section do we belong to
    2910             :  *                   1: Standard monadic or dyadic OPERATORS
    2911             :  *                   2: CONVERSIONS and similar elementary functions
    2912             :  *                   3: functions related to COMBINATORICS
    2913             :  *                   4: TRANSCENDENTAL functions, etc.
    2914             :  *   char *code   : GP prototype, aka Parser Code (see libpari's manual)
    2915             :  *                  if NULL, use valence instead.
    2916             :  *   char *help   : short help text (init to NULL).
    2917             :  *   void *pvalue : push_val history.
    2918             :  *   long arity   : maximum number of arguments.
    2919             :  *   entree *next : next entree (init to NULL, used in hashing code). */
    2920             : #include "init.h"
    2921             : #include "default.h"

Generated by: LCOV version 1.13