1 /* -*- mode:C++ ; compile-command: "g++-3.4 -I.. -g -c pari.cc" -*- */
2 #include "giacPCH.h"
3 
4 /*  PARI interface
5  *  Copyright (C) 2001,14 B. Parisse, Institut Fourier, 38402 St Martin d'Heres
6  *
7  *  This program is free software; you can redistribute it and/or modify
8  *  it under the terms of the GNU General Public License as published by
9  *  the Free Software Foundation; either version 3 of the License, or
10  *  (at your option) any later version.
11  *
12  *  This program is distributed in the hope that it will be useful,
13  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  *  GNU General Public License for more details.
16  *
17  *  You should have received a copy of the GNU General Public License
18  *  along with this program. If not, see <http://www.gnu.org/licenses/>.
19  */
20 #if !defined NSPIRE
21 using namespace std;
22 #endif
23 #ifdef HAVE_CONFIG_H
24 #include "config.h"
25 #endif
26 
27 #include <stdexcept>
28 #include "gen.h"
29 #include "identificateur.h"
30 #include "sym2poly.h"
31 #include "plot.h"
32 #include "prog.h"
33 #include "usual.h"
34 #include "input_lexer.h"
35 #include "modpoly.h"
36 #include "giacintl.h"
37 #ifdef USE_GMP_REPLACEMENTS
38 #undef HAVE_LIBPARI
39 #endif
40 
41 #ifdef HAVE_LIBPARI
42 #ifdef HAVE_PTHREAD_H
43 #include <pthread.h>
44 #endif
45 
abs(long int & l)46 static long int abs(long int & l){
47   if (l<0)
48     return -l;
49   else
50     return l;
51 }
52 #include "pari.h"
53 extern "C" {
54 #include <pari/pari.h>
55 #include <pari/paripriv.h>
56 #ifdef ENABLE_TLS
57   extern THREAD void *PARI_stack_limit;
58 #else
59   extern void *PARI_stack_limit;
60 #endif
61   extern entree functions_basic[];
62 }
63 #if PARI_VERSION_CODE<PARI_VERSION(2,4,0) // 132096
64 #define PARI23
65 #else
66 jmp_buf env;
67 static void
gp_err_recover(long numerr)68 gp_err_recover(long numerr)
69 {
70   longjmp(env, numerr);
71 }
72 
73 #endif
74 #include <cstdlib>
75 
76 #ifndef NO_NAMESPACE_GIAC
77 namespace giac {
78 #endif // ndef NO_NAMESPACE_GIAC
79 
80 #ifdef HAVE_LIBPTHREAD
81   static pthread_mutex_t * pari_mutex_ptr = 0;
82 #endif
83 
84 #if PARI_VERSION_CODE >= PARI_VERSION(2,8,0)
85 #define PARI_DYNAMIC_STACK
86 #endif
87 
88   static map<string,entree *> pari_function_table;
do_giac_pari_init(long maxprime)89   static void do_giac_pari_init(long maxprime){
90 #ifdef PARI_DYNAMIC_STACK
91     long pari_mem_size=512000;
92 #else
93     long pari_mem_size=64000000;
94 #endif
95     if (getenv("PARI_SIZE")){
96       string pari_size_s(getenv("PARI_SIZE"));
97       pari_mem_size= atoi(pari_size_s.c_str());
98     }
99     // do not initialize INIT_JMP so that PARI error do not exit
100     pari_init_opts(pari_mem_size,maxprime,INIT_SIGm | INIT_DFTm);
101 #ifdef PARI_DYNAMIC_STACK
102     paristack_setsize(pari_mem_size, (1<<30)); // pari 2.8
103     //paristack_alloc(pari_mem_size, (1<<30)); // pari 2.8
104 #endif
105     entree * ptr=functions_basic;
106     for (;ptr->name;++ptr){
107       pari_function_table[ptr->name]=ptr;
108     }
109     // initialize variable ordering x,y,z
110     gp_read_str("[x,y,z,t]");
111   }
112 
113 #if 0
114   struct giac_pari_init {
115     giac_pari_init(long maxprime) {
116       if(!avma){
117 	do_giac_pari_init(maxprime);
118       }
119     }
120   };
121   static long pari_maxprime=100000;
122   static giac_pari_init bidon(pari_maxprime);
123 #else
124   static long pari_maxprime=100000;
get_pari_avma()125   long get_pari_avma() {
126     if(!avma){
127       do_giac_pari_init(pari_maxprime);
128     }
129     return avma;
130   }
131 #endif
132 
133   static gen pow2sizeof_long(pow(256,sizeof(long)));
134   // Conversion of a GEN integer to a gen, using Horner method
t_INT2gen(const GEN & G)135   static gen t_INT2gen(const GEN & G){
136     long Gs=signe(G);
137     if (!Gs)
138       return 0;
139     //setsigne(G,1);
140     long Gpl=lgefint(G)-2;
141     // use one of the following code depending how pari codes long integers
142     // Code 1
143 #if !defined(__APPLE__) && !defined(WIN32) || defined(WIN64)
144     ref_mpz_t * mz = new ref_mpz_t;
145     mpz_realloc2(mz->z,32*Gpl);
146     mpz_import(mz->z,Gpl,-1,sizeof(GEN),0,0,&G[2]);
147     //setsigne(G,Gs);
148     if (Gs>0)
149       return mz;
150     else
151       return -gen(mz);
152 #else
153     // Code 2
154     --Gpl;
155     long * Gp=int_MSW(G);
156     gen res;
157     for (int i=0;i<=Gpl;++i){
158 #ifdef INT128
159       res=res*pow2sizeof_long+int128_t((ulonglong)*Gp);
160 #else
161       res=res*pow2sizeof_long+longlong(unsigned(*Gp));
162 #endif
163       Gp=int_precW(Gp);
164     }
165     return Gs<0?-res:res;
166 #endif
167   }
168 
t_REAL2gen(const GEN & G)169   static gen t_REAL2gen(const GEN & G){
170     long Gs=signe(G);
171     if (!Gs)
172       return 0.0;
173     long n=lg(G);
174     gen res;
175     for (int i=2;i<n;++i){
176 #ifdef x86_64 // FIXME make a gen constructor from ulonglong
177       unsigned long u=G[i];
178       res=res*pow2sizeof_long;
179       if (u%2)
180 	res += 1;
181       res += 2*gen(longlong(u>>1));
182 #else
183       res=res*pow2sizeof_long+longlong(unsigned(G[i]));
184 #endif
185     }
186     res=res*pow(plus_two,int(expo(G)+1-bit_accuracy(n)));
187     res=_evalf(makesequence(res,int(n/3.3)),context0);
188     return Gs<0?-res:res;
189   }
190 
t_POL2gen(const GEN & G,const vecteur & vars)191   static gen t_POL2gen(const GEN & G,const vecteur & vars){
192     if (!signe(G))
193       return 0;
194     long n=lg(G);
195     vecteur res;
196     for (long i=2;i<n;++i){
197       res.push_back(GEN2gen((GEN)G[i],vars));
198     }
199     reverse(res.begin(),res.end());
200     long vn=varn(G);
201     gen x;
202     if (vn<long(vars.size())){
203       x=vars[vn];
204       return symb_horner(res,x);
205     }
206     if (!vars.empty() && vn==long(vars.size())){
207       x=vars[vn-1];
208       return symb_horner(res,x);
209     }
210     else
211       return gen(res,_POLY1__VECT);
212   }
213 
t_VEC2gen(const GEN & G,const vecteur & vars,long debut,long fin)214   static gen t_VEC2gen(const GEN & G,const vecteur & vars,long debut,long fin){
215     vecteur res;
216     for (long i=debut;i<fin;++i){
217       res.push_back(GEN2gen((GEN)G[i],vars));
218     }
219     return gen(res,typ(G)==t_VEC?0:99);
220   }
221 
t_VECSMALL2gen(const GEN & G)222   static gen t_VECSMALL2gen(const GEN & G){
223     long fin=lg(G);
224     vecteur res;
225     for (long i=1;i<fin;++i){
226       res.push_back((int)G[i]);
227     }
228     return gen(res,98);
229   }
230 
t_MOD2gen(const GEN & G,const vecteur & vars)231   static gen t_MOD2gen(const GEN & G,const vecteur & vars){
232     return makemod(GEN2gen((GEN) G[2],vars),GEN2gen((GEN) G[1],vars));
233   }
234 
t_POLMOD2gen(const GEN & G,const vecteur & vars)235   static gen t_POLMOD2gen(const GEN & G,const vecteur & vars){
236     gen tmp=_fxnd(GEN2gen((GEN) G[2],vars),context0);
237     gen n(tmp),d(1);
238     if (tmp.type==_VECT && tmp._VECTptr->size()==2){
239       n=tmp._VECTptr->front();
240       if (n.type==_VECT) n.subtype=_POLY1__VECT;
241       d=tmp._VECTptr->back();
242     }
243     return eval(symbolic(at_rootof,makesequence(n,change_subtype(GEN2gen((GEN) G[1],vars),_POLY1__VECT))),1,context0)/d;
244     find_or_make_symbol("Mod",tmp,0,false,context0);
245     return symbolic(at_of,makesequence(tmp,gen(makevecteur(GEN2gen((GEN) G[2],vars),GEN2gen((GEN) G[1],vars)),_SEQ__VECT)));
246   }
247 
t_COMPLEX2gen(const GEN & G,const vecteur & vars)248   static gen t_COMPLEX2gen(const GEN & G,const vecteur & vars){
249     return gen(GEN2gen((GEN) G[1],vars),GEN2gen((GEN) G[2],vars));
250   }
251 
t_FRAC2gen(const GEN & G,const vecteur & vars)252   static gen t_FRAC2gen(const GEN & G,const vecteur & vars){
253     return fraction(GEN2gen((GEN)G[1],vars),GEN2gen((GEN)G[2],vars));
254   }
255 
t_QUAD2gen(const GEN & G,const vecteur & vars)256   static gen t_QUAD2gen(const GEN & G,const vecteur & vars){
257     // use w__IDNT_e like pari for all quadratics
258     return GEN2gen((GEN) G[1],vars)+w__IDNT_e*GEN2gen((GEN)G[2],vars);
259   }
260 
t_PADIC2gen(const GEN & G,const vecteur & vars)261   static gen t_PADIC2gen(const GEN & G,const vecteur & vars){
262     gen O;
263     find_or_make_symbol("O",O,0,false,context0);
264     gen p(GEN2gen((GEN) G[2],vars)),val(longlong(valp(G)));
265     return pow(p,val,context0)*(GEN2gen((GEN) G[4],vars)+symbolic(at_of,makesequence(O,symb_pow(p,longlong(precp(G)))))); // removed symb_quote for p-adic
266   }
267 
268   // WARNING: If g is a matrix this print the transpose of the matrix
GEN2string(const GEN & g)269   static string GEN2string(const GEN & g){
270     // cerr << typ(g) << " " << t_MAT << '\n';
271     char * ch;
272     string s;
273     if ((typ(g)==t_MAT) || (typ(g)==t_COL)){
274       int taille=lg(g);
275       s +="[";
276       for (int i=1;i<taille;++i){
277 	s += GEN2string((long *)g[i]);
278 	if (i==taille-1)
279 	  s+="]";
280 	else
281 	  s+=",";
282       }
283       return s;
284     }
285     ch=GENtostr(g);
286     s=ch;
287 #ifndef HAVE_LIBGC
288     free(ch);
289 #endif
290     return s;
291   }
292 
default2gen(const GEN & G)293   static gen default2gen(const GEN &G){
294     string s=GEN2string(G);
295     gen g;
296     try {
297       g=gen(s,context0);
298     } catch(...){
299       return string2gen(s,false);
300     }
301     return g;
302   }
303 
GEN2gen(const GEN & G,const vecteur & vars)304   gen GEN2gen(const GEN & G,const vecteur & vars){
305 #ifdef EMCC
306     return default2gen(G);
307 #endif
308     switch (typ(G)){
309     case t_INT:
310       return t_INT2gen(G);
311     case t_INTMOD:
312       return t_MOD2gen(G,vars);
313     case t_POLMOD:
314       return t_POLMOD2gen(G,vars);
315     case t_FRAC: case t_RFRAC:
316       return t_FRAC2gen(G,vars);
317     case t_COMPLEX:
318       return t_COMPLEX2gen(G,vars);
319     case t_REAL:
320       return t_REAL2gen(G);
321     case t_POL:
322       return t_POL2gen(G,vars);
323     case t_VEC: case t_COL:
324       return t_VEC2gen(G,vars,1,lg(G));
325     case t_VECSMALL:
326       return t_VECSMALL2gen(G);
327     case t_MAT:
328       return _tran(t_VEC2gen(G,vars,1,lg(G)),context0);
329     case t_LIST:
330 #ifdef PARI23
331       return t_VEC2gen(G,vars,2,lgeflist(G));
332 #else
333       return t_VEC2gen(G,vars,2,list_nmax(G));
334 #endif
335     case t_STR:
336       return string2gen(GSTR(G),false);
337     case t_QUAD:
338       return t_QUAD2gen(G,vars);
339     case t_PADIC:
340       return t_PADIC2gen(G,vars);
341     default:
342       return default2gen(G);
343     }
344   }
345 
346   static std::string pariprint(const gen & e,int varnum,GIAC_CONTEXT);
347 
pariprint_VECT(const vecteur & v,int varnum,int subtype,GIAC_CONTEXT)348   static string pariprint_VECT(const vecteur & v,int varnum,int subtype,GIAC_CONTEXT){
349     string s;
350     const_iterateur it=v.begin(),itend=v.end();
351     if (subtype==_POLY1__VECT){
352       if (v.empty())
353 	return "0";
354       string tmp=")*x"+print_INT_(varnum)+"+";
355       ++varnum;
356       s=pariprint(*it,varnum,contextptr);
357       for (++it;it!=itend;++it){
358 	s="("+s+tmp+pariprint(*it,varnum,contextptr);
359       }
360       --varnum;
361       return s;
362     }
363     if (subtype!=_SEQ__VECT)
364       s="[";
365     for (;it!=itend;++it){
366       s += pariprint(*it,varnum,contextptr);
367       if (it+1!=itend)
368 	s += ",";
369     }
370     if (subtype!=_SEQ__VECT)
371       s+="]";
372     return s;
373   }
374 
pariprintmatrice(const gen & e,int varnum,GIAC_CONTEXT)375   static string pariprintmatrice(const gen & e,int varnum,GIAC_CONTEXT){
376     string res = "[";
377     const_iterateur it=e._VECTptr->begin(),itend=e._VECTptr->end();
378     for (;it!=itend;++it){
379       res += pariprint_VECT(*it->_VECTptr,varnum,_SEQ__VECT,contextptr);
380       if (it+1!=itend)
381 	res += ";";
382       else
383 	res +="]";
384     }
385     return res;
386   }
387 
pariprint(const gen & e,int varnum,GIAC_CONTEXT)388   static string pariprint(const gen & e,int varnum,GIAC_CONTEXT){
389     int py=python_compat(contextptr);
390     python_compat(0,contextptr);
391     int save_maple_mode=xcas_mode(contextptr);
392     xcas_mode(contextptr)=1;
393     string res;
394     switch (e.type){
395     case _INT_:
396       res=print_INT_(e.val);
397       break;
398     case _ZINT:
399       res=e.print(contextptr);
400       break;
401     case _CPLX:
402       res=e._CPLXptr->print(contextptr)+"+I*"+(e._CPLXptr+1)->print(contextptr);
403       break;
404     case _VECT:
405       if (ckmatrix(e))
406 	res=pariprintmatrice(e,varnum,contextptr);
407       else
408 	res=pariprint_VECT(*e._VECTptr,varnum,e.subtype,contextptr);
409       break;
410     case _SYMB:
411       if (e._SYMBptr->sommet==at_rootof && e._SYMBptr->feuille.type==_VECT && e._SYMBptr->feuille._VECTptr->size()==2)
412 	res="Mod("+pariprint(e._SYMBptr->feuille._VECTptr->front(),varnum,contextptr)+","+pariprint(e._SYMBptr->feuille._VECTptr->back(),varnum,contextptr)+")";
413       else {
414 	if (e._SYMBptr->sommet.ptr()->printsommet) // FIXME
415 	  res=e.print(contextptr);
416 	else
417 	  res=string(e._SYMBptr->sommet.ptr()->s)+"("+pariprint(e._SYMBptr->feuille,varnum,contextptr)+")";
418       }
419       break;
420     case _MOD:
421       res= "Mod("+pariprint(*e._MODptr,varnum,contextptr)+","+pariprint(*(e._MODptr+1),varnum,contextptr)+")";
422       break;
423     case _FRAC:
424       res= pariprint(*e._FRACptr,varnum,contextptr)+"/("+pariprint(*(e._FRACptr+1),varnum,contextptr)+")";
425       break;
426     default: // _SYMB with printsommetasoperator, _DOUBLE_, _REAL, _IDNT
427       res=e.print(contextptr);
428     }
429     xcas_mode(contextptr)=save_maple_mode;
430     python_compat(py,contextptr);
431     return res;
432   }
433 
zint2GEN(const gen & g)434   static GEN zint2GEN(const gen & g){
435     mpz_t * zz=g._ZINTptr;
436     int sgn=mpz_sgn(*zz);
437     if (!sgn)
438       return utoi(0);
439     int count=mpz_sizeinbase(*zz,2);
440     if (count % (8*sizeof(GEN)))
441       count=count/(8*sizeof(GEN))+3;
442     else
443       count=count/(8*sizeof(GEN))+2;
444     GEN G=cgetg(count,t_INT);
445     //size_t countp;
446     // mpz_export(&G[2],&countp,-1,sizeof(GEN),0,0,*zz);
447     setlgefint(G,count);
448     setsigne(G,sgn);
449     // return G;
450     gen q(abs(g)),tmp,r;
451     GEN Gstep(int2n(16));
452     vector<int> v;
453     for (;!is_zero(q);){
454       r=irem(q,gen(1<<16),tmp);
455       v.push_back(r.val);
456       q=tmp;
457     }
458     int s=v.size();
459     GEN res(utoi(0));
460     for (int i=s-1;i>=0;--i){
461       res=gmul(res,Gstep);
462       res=gadd(res,utoi(v[i]));
463     }
464     setsigne(res,sgn);
465     return res;
466   }
467 
468   static GEN ingen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT);
469 
vect2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT)470   static GEN vect2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT){
471     vecteur v (*g._VECTptr);
472     int n=v.size(),decal=1;
473     GEN res;
474     if (g.subtype==_POLY1__VECT){
475       decal=2;
476       res=cgetg(n+decal,t_POL);
477       reverse(v.begin(),v.end());
478     }
479     else {
480       res=cgetg(n+decal,g.subtype==99?t_COL:(g.subtype==98?t_VECSMALL:t_VEC));
481     }
482     for (int i=0;i<n;++i)
483       gel(res,i+decal)=ingen2GEN(v[i],vars,contextptr);
484     if (decal==2){
485       setsigne(res,1);
486       setvarn(res,0);
487     }
488     return res;
489   }
490 
mat2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT)491   static GEN mat2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT){
492     matrice M = mtran(*g._VECTptr);
493     int n=M.size(),m=M[0]._VECTptr->size();
494     GEN res=cgetg(n+1,t_MAT);
495     for (int i=1;i<=n;++i){
496       GEN resi=gel(res,i)=cgetg(m+1,t_COL);
497       vecteur & v = *M[i-1]._VECTptr;
498       for (int j=1;j<=m;++j){
499 	gel(resi,j)=ingen2GEN(v[j-1],vars,contextptr);
500       }
501     }
502     return res;
503   }
504 
cplx2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT)505   static GEN cplx2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT){
506     GEN res=cgetg(3,t_COMPLEX);
507     gel(res,1)=ingen2GEN(*g._CPLXptr,vars,contextptr);
508     gel(res,2)=ingen2GEN(*(g._CPLXptr+1),vars,contextptr);
509     return res;
510   }
511 
frac2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT)512   static GEN frac2GEN(const gen & g,const vecteur & vars,GIAC_CONTEXT){
513     GEN res=cgetg(3,t_FRAC);
514     gel(res,1)=ingen2GEN(g._FRACptr->num,vars,contextptr);
515     gel(res,2)=ingen2GEN(g._FRACptr->den,vars,contextptr);
516     return res;
517   }
518 
real2GEN(const gen & e_,GIAC_CONTEXT)519   static GEN real2GEN(const gen & e_,GIAC_CONTEXT){
520     gen e(e_);
521     bool neg=false;
522     if (is_strictly_positive(-e,contextptr)){
523       e=-e;
524       neg=true;
525     }
526     int prec=53;
527 #ifdef HAVE_LIBMPFR
528     if (e.type==_REAL)
529       prec=mpfr_get_prec(e._REALptr->inf);
530 #endif
531     string s(e.print(contextptr));
532     GEN g=strtor(s.c_str(),prec);
533     if (neg)
534       g=gneg(g);
535     if (debug_infolevel)
536       CERR << "real converted to pari " << GEN2gen(g,vecteur(0)) << '\n';
537     else e=GEN2gen(g,vecteur(0));
538     // for some strange reason, converting g to a gen fixes a bug in conversion
539     return g;
540   }
541 
ingen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT)542   static GEN ingen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT){
543     switch (e.type){
544     case _INT_:
545       return stoi(e.val);
546     case _ZINT:
547       return zint2GEN(e);
548     case _CPLX:
549       return cplx2GEN(e,vars,contextptr);
550     case _FRAC:
551       return frac2GEN(e,vars,contextptr);
552     case _DOUBLE_: case _REAL:
553       return real2GEN(e,contextptr);
554     case _VECT:
555       if (ckmatrix(e))
556 	return mat2GEN(e,vars,contextptr);
557       else
558 	return vect2GEN(e,vars,contextptr);
559     }
560     // add vars to e
561     string s=pariprint(e,0,contextptr);
562     vecteur vars_(vars);
563     // if (vars_.empty()) vars_.push_back(vx_var);
564     if (!vars_.empty())
565       s="["+(vars_.size()==1?vars_.front().print():print_VECT(vars_,_SEQ__VECT,contextptr))+","+s+"]";
566     GEN res= gp_read_str((char *) s.c_str());
567     // gp_read_str seems to have problems with large strings (s.size()>2^16?)
568     return vars_.empty()?res:gel(res,1+vars_.size());
569   }
gen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT)570   GEN gen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT){
571 #ifdef PARI23
572     if (setjmp(GP_DATA->env)){
573       setsizeerr(gettext("Error in PARI subsystem"));
574     }
575 #else
576     cb_pari_err_recover=gp_err_recover;
577     if (setjmp(env)){
578       setsizeerr(gettext("Error in PARI subsystem"));
579     }
580 #endif
581     return ingen2GEN(e,vars,contextptr);
582   }
583 
pari_cleanup(void * arg)584   static void pari_cleanup(void * arg) {
585 #ifdef HAVE_LIBPTHREAD
586     if (arg)
587       pthread_mutex_unlock((pthread_mutex_t *)arg);
588 #endif
589   }
590 
check_pari_mutex()591   int check_pari_mutex(){
592 #ifdef HAVE_LIBPTHREAD
593     if (!pari_mutex_ptr){
594       pthread_mutex_t tmp=PTHREAD_MUTEX_INITIALIZER;
595       pari_mutex_ptr=new pthread_mutex_t(tmp);
596     }
597     return pthread_mutex_trylock(pari_mutex_ptr);
598 #else
599     return 0;
600 #endif
601   }
602 
abort_if_locked()603   void abort_if_locked(){
604     if (check_pari_mutex())
605       setsizeerr(gettext("PARI locked by another thread. Try again later.\nIf PARI is locked by a cancelled thread, you can unlock it with pari_unlock()"));
606   }
607 
pari_isprime(const gen & e,int certif)608   gen pari_isprime(const gen & e,int certif){
609     gen tmp;
610     abort_if_locked();
611 #ifdef HAVE_LIBPTHREAD
612     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
613 #endif
614     long av=get_pari_avma();
615     tmp=GEN2gen(gisprime(gen2GEN(e,vecteur(0),0),certif),vecteur(0));
616     avma=av;
617 #ifdef HAVE_LIBPTHREAD
618     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
619     pthread_cleanup_pop(0);
620 #endif
621     return tmp;
622   }
623 
pari_ifactor(const gen & e)624   string pari_ifactor(const gen & e){
625     abort_if_locked();
626     string s;
627 #ifdef HAVE_LIBPTHREAD
628     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
629 #endif
630     long av=get_pari_avma();
631     GEN g=gen2GEN(e,vecteur(0),0);
632     GEN gf=factorint(g,0);
633     s=GEN2string(gf);
634     avma=av;
635 #ifdef HAVE_LIBPTHREAD
636     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
637     pthread_cleanup_pop(0);
638 #endif
639     return s;
640   }
641 
pari_gamma(const gen & e)642   gen pari_gamma(const gen & e){
643     abort_if_locked();
644     gen res;
645 #ifdef HAVE_LIBPTHREAD
646     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
647 #endif
648     long av=get_pari_avma();
649     GEN g=gen2GEN(e,vecteur(0),0);
650     GEN gf=ggamma(g,precision(g));
651     res=GEN2gen(gf,vecteur(0));
652     avma=av;
653 #ifdef HAVE_LIBPTHREAD
654     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
655     pthread_cleanup_pop(0);
656 #endif
657     return res;
658   }
659 
pari_zeta(const gen & e)660   gen pari_zeta(const gen & e){
661     abort_if_locked();
662     gen res;
663 #ifdef HAVE_LIBPTHREAD
664     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
665 #endif
666     long av=get_pari_avma();
667     GEN g=gen2GEN(e,vecteur(0),0);
668     GEN gf=gzeta(g,precision(g));
669     res=GEN2gen(gf,vecteur(0));
670     avma=av;
671 #ifdef HAVE_LIBPTHREAD
672     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
673     pthread_cleanup_pop(0);
674 #endif
675     return res;
676   }
677 
pari_psi(const gen & e)678   gen pari_psi(const gen & e){
679     abort_if_locked();
680     gen res;
681 #ifdef HAVE_LIBPTHREAD
682     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
683 #endif
684     long av=get_pari_avma();
685     GEN g=gen2GEN(e,vecteur(0),0);
686     GEN gf=gpsi(g,precision(g));
687     res=GEN2gen(gf,vecteur(0));
688     avma=av;
689 #ifdef HAVE_LIBPTHREAD
690     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
691     pthread_cleanup_pop(0);
692 #endif
693     return res;
694   }
695 
pari_ffinit(const gen & p,int n)696   gen pari_ffinit(const gen & p,int n){
697     gen tmp;
698 #ifdef HAVE_LIBPTHREAD
699     abort_if_locked();
700     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
701 #endif
702     long av=get_pari_avma();
703     tmp=GEN2gen(ffinit(gen2GEN(p,vecteur(0),0),n,0),vecteur(0));
704     avma=av;
705 #ifdef HAVE_LIBPTHREAD
706     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
707     pthread_cleanup_pop(0);
708 #endif
709     return tmp;
710   }
711 
712   // for factorization over Z when many modular factors arise
713   // This is a call to PARI combine_factors
714   // WARNING: You must remove static from the declaration of combine_factors
715   // in pari/src/basemath/polarith2.c
716   // GEN combine_factors(GEN a, GEN famod, GEN p, long klim, long hint);
pari_lift_combine(const vecteur & a,const vector<vecteur> & factmod,gen & modulo,vector<vecteur> & res)717   bool pari_lift_combine(const vecteur & a,const vector<vecteur> & factmod,gen & modulo,vector<vecteur> & res){
718 #ifdef PARI23
719     long av=get_pari_avma();
720     GEN pari_a=gen2GEN(r2e(a,x__IDNT_e,context0),vecteur(0),0);
721     string s("[");
722     vector<vecteur>::const_iterator it=factmod.begin(),itend=factmod.end();
723     for (;it!=itend;){
724       s += r2e(*it,x__IDNT_e,context0).print();
725       ++it;
726       if (it==itend)
727 	break;
728       s+=",";
729     }
730     s+="]";
731     // cerr << s << '\n';
732     GEN pari_factmod=gp_read_str((char *) s.c_str());
733     GEN pari_modulo=gen2GEN(modulo,vecteur(0),0);
734     GEN pari_res=combine_factors(pari_a,pari_factmod,pari_modulo,0,1);
735     // back conversion
736     string res_s=GENtostr(pari_res);
737     gen res_v(res_s.substr(0,res_s.size()-1),context0);
738     if (res_v.type!=_VECT)
739       setsizeerr();
740     const_iterateur jt=res_v._VECTptr->begin(),jtend=res_v._VECTptr->end();
741     for (;jt!=jtend;++jt){
742       res.push_back(*e2r(*jt,x__IDNT_e,context0)._VECTptr);
743     }
744     avma=av;
745     return true;
746 #else
747     return false;
748 #endif
749   }
750 
pari_exec(const string & s,GIAC_CONTEXT)751   static gen pari_exec(const string & s,GIAC_CONTEXT){
752     long av=get_pari_avma();
753     void * save_pari_stack_limit = PARI_stack_limit;
754     PARI_stack_limit=0; // required since the stack changed
755 #ifdef PARI23
756     if (setjmp(GP_DATA->env))
757 #else
758     cb_pari_err_recover=gp_err_recover;
759     if (setjmp(env))
760 #endif
761       {
762 #ifdef HAVE_LIBPTHREAD
763 	if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
764 #endif
765 	avma = av;
766 	*logptr(contextptr) << gettext("Error in PARI subsystem") << '\n';
767 	PARI_stack_limit = save_pari_stack_limit ;
768 	// setsizeerr();
769 	return undef;
770       }
771     GEN gres= gp_read_str((char *) s.c_str());
772     gen res=GEN2gen(gres,vecteur(0));
773     avma=av;
774     PARI_stack_limit = save_pari_stack_limit ;
775     return res;
776   }
777 
778 #include "input_parser.h"
779 #define _ARGS_ argvec[0], argvec[1], argvec[2], argvec[3],\
780                argvec[4], argvec[5], argvec[6], argvec[7], argvec[8]
781   // args=pari_function_name, arg1, ...
782   // or pari_function_name quoted to define a function
783 
784   enum {
785     RET_GEN=0,
786     RET_VOID=1,
787     RET_INT=2,
788     RET_LONG=3
789   };
790   typedef GEN (*PFGEN)(ANYARG);
791 
792   extern const unary_function_ptr * const  at_pari;
in_pari(const gen & args,GIAC_CONTEXT)793   static gen in_pari(const gen & args,GIAC_CONTEXT){
794     if (args.type<_IDNT)
795       return args;
796     vecteur v(gen2vecteur(args));
797     int vs=v.size();
798     if (!vs){ // export all pari functions
799       entree * ptr=functions_basic;
800       gen tmp; int lextype;
801       string redef;
802       for (;ptr->name;++ptr){
803 	pari_function_table[ptr->name]=ptr;
804 	lextype=find_or_make_symbol(ptr->name,tmp,0,false,contextptr);
805 	if (lextype==T_SYMBOL)
806 	  sto(symbolic(at_pari,string2gen(ptr->name,false)),tmp,contextptr);
807 	else
808 	  redef += string(ptr->name) + " ";
809 	find_or_make_symbol(string("pari_")+ptr->name,tmp,0,false,contextptr);
810 	sto(symbolic(at_pari,string2gen(ptr->name,false)),tmp,contextptr);
811       }
812       return string2gen("All PARI functions are now defined with the pari_ prefix.\nPARI functions are also defined without prefix except:\n"+redef+"\nWhen working with p-adic numbers use them in a pari() call\nType ?pari for short help\nInside xcas, try Help->Manuals->PARI for HTML help",false);
813     }
814     if (v[0].is_symb_of_sommet(at_quote)){
815       if (vs==1)
816 	return symbolic(at_pari,args);
817       v[0]=v[0]._SYMBptr->feuille;
818     }
819     for (int i=1;i<vs;i++)
820       v[i]=v[i].eval(eval_level(contextptr),contextptr);
821     vecteur vars(1,identificateur("O"));
822     lidnt(v,vars,false);
823     vars.erase(vars.begin());
824     bool parse_all=false;
825     long av=get_pari_avma();
826 #ifdef PARI23
827     if (setjmp(GP_DATA->env)){
828       avma = av;
829       parse_all=true;
830       if (setjmp(GP_DATA->env)){ // if (setjmp(GP_DATA->env)){
831 #ifdef HAVE_LIBPTHREAD
832 	if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
833 #endif
834 	avma = av;
835 	*logptr(contextptr) << gettext("Error in PARI subsystem") << '\n';
836 	// setsizeerr();
837 	return undef;
838       }
839     }
840 #else
841     cb_pari_err_recover=gp_err_recover;
842     if (setjmp(env)){
843       avma = av;
844       parse_all=true;
845       if (setjmp(env)){ // if (setjmp(GP_DATA->env)){
846 #ifdef HAVE_LIBPTHREAD
847 	if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
848 #endif
849 	avma = av;
850 	*logptr(contextptr) << gettext("Error in PARI subsystem") << '\n';
851 	// setsizeerr();
852 	return undef;
853       }
854     }
855 #endif
856 #ifdef EMCC
857     parse_all=true;
858     // otherwise fibonacci or lngamma do not work, probably because of the 'L'
859     // conversion (long argument expected) or precision (same)
860 #endif
861     if (!parse_all && v[0].type==_STRNG) {
862       string vstr=*v[0]._STRNGptr;
863       if (vstr!="") {
864 	void * save_pari_stack_limit = PARI_stack_limit;
865 	PARI_stack_limit=0; // required since the stack changed
866 	if (vs==1)
867 	  return symbolic(at_pari,args);
868 	map<string,entree *>::const_iterator i = pari_function_table.find(vstr);
869 	// look at function prototype for return value
870 	// and call code, from anal.c line around 1990
871 	unsigned int ret;
872 	if (i!=pari_function_table.end()){
873 	  const char * s =i->second->code;
874 	  if      (*s <  'a')   ret = RET_GEN;
875 	  else if (*s == 'v') { ret = RET_VOID; s++; }
876 	  else if (*s == 'i') { ret = RET_INT;  s++; }
877 	  else if (*s == 'l') { ret = RET_LONG; s++; }
878 	  else                  ret = RET_GEN;
879 	  void * call= i->second->value;
880 	  // translate gen to GEN's
881 	  GEN argvec[9]={0,0,0,0,0,0,0,0,0},res; long m;
882 	  int k=0;
883 	  for (int j=1;k<9 && *s && *s!='\n';++s){
884 	    switch(*s){
885 	    case 'L': // long
886 	      if (j==vs) {
887 		gen res=string2gen("PARI: Bad argument count",false);
888 		res.subtype=-1;
889 		return res;
890 	      }
891 	      argvec[k]= (GEN) v[j].val;
892 	      ++j; ++k;
893 	      break;
894 	    case 'P': // default precision
895 	      argvec[k] = (GEN) precdl; k++; break;
896 	    case 'D': //default param
897 	      {
898 		++s;
899 		switch(*s){
900 		case 'G': case '&': case 'I': case 'V':
901 		  if (j<vs)
902 		    argvec[k]=ingen2GEN(v[j],vars,contextptr);
903 #ifdef PARI23
904 		  else
905 		    argvec[k]=utoi(0);
906 #endif
907 		  ++j; ++k;
908 		  break;
909 		case 'n':
910 		  if (j<vs){
911 		    int pos=equalposcomp(vars,v[j]);
912 		    if (pos)
913 		      argvec[k]=(long int*)(pos -1);
914 		  }
915 		  else
916 		    argvec[k]=0;
917 		  ++j; ++k;
918 		  break;
919 		default:
920 		  if (j<vs)
921 		    argvec[k]=(long int*) v[j].val;
922 		  else
923 		    argvec[k]=0;
924 		  ++k; ++j;
925 		  while (*s!= ',') s++;
926 		  s++;
927 		  while (*s!= ',') s++;
928 		}
929 		break;
930 	      }
931 	    case 'p':
932 	      argvec[k]=(GEN) precreal;
933 	      ++k;
934 	      break;
935 	    default:
936 	      if (j==vs) {
937 		gen res=string2gen("PARI: Bad argument count",false);
938 		res.subtype=-1;
939 		return res;
940 	      }
941 	      argvec[k]=ingen2GEN(v[j],vars,contextptr);
942 	      ++j; ++k;
943 	      break;
944 	    }
945 	  }
946 	  switch (ret)
947 	    {
948 	    case RET_GEN:
949 	      res = ((PFGEN)call)(_ARGS_);
950 	      break;
951 
952 	    case RET_INT:
953 	      m = (long)((int (*)(ANYARG))call)(_ARGS_);
954 	      res = stoi(m); break;
955 
956 	    case RET_LONG:
957 	      m = ((long (*)(ANYARG))call)(_ARGS_);
958 	      res = stoi(m); break;
959 
960 	    case RET_VOID:
961 	      ((void (*)(ANYARG))call)(_ARGS_);
962 	      res = gnil; break;
963 	    }
964 	  // cerr << GEN2string(res) << '\n';
965 	  gen resg(GEN2gen(res,vars));
966 	  PARI_stack_limit = save_pari_stack_limit ;
967 	  avma=av;
968 	  return resg;
969 	} // end if (i!=pari_function_table.end())
970       } // end if vstr!=""
971       if (vstr=="" && vs==2){
972 	long av=get_pari_avma();
973 	gen res= GEN2gen(gen2GEN(v[1],vars,contextptr),vars);
974 	avma=av;
975 	return res;
976       }
977     } // end if (!parse_all ...)
978     string s;
979     if (v[0].type==_FUNC)
980       s=v[0]._FUNCptr->ptr()->s;
981     else
982       s=gen2string(v[0]);
983     if (vs>1){
984       s+="(";
985       for (int i=1;i<vs;){
986 	s += pariprint(v[i],0,contextptr);
987 	++i;
988 	if (i==vs)
989 	  break;
990 	s += ",";
991       }
992       s +=")";
993     }
994     return pari_exec(s,contextptr);
995   }
_pari(const gen & args,GIAC_CONTEXT)996   gen _pari(const gen & args,GIAC_CONTEXT){
997     if ( args.type==_STRNG && args.subtype==-1) return  args;
998     abort_if_locked();
999     gen res;
1000 #ifdef HAVE_LIBPTHREAD
1001     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
1002 #endif
1003     res=in_pari(args,contextptr);
1004 #ifdef HAVE_LIBPTHREAD
1005     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
1006     pthread_cleanup_pop(0);
1007 #endif
1008     return res;
1009   }
1010   static const char _pari_s []="pari";
1011   static define_unary_function_eval (__pari,&_pari,_pari_s);
1012   define_unary_function_ptr5( at_pari ,alias_at_pari,&__pari,_QUOTE_ARGUMENTS,true);
1013 
_pari_unlock(const gen & args,GIAC_CONTEXT)1014   gen _pari_unlock(const gen & args,GIAC_CONTEXT){
1015     if ( args.type==_STRNG && args.subtype==-1) return  args;
1016     int locked=check_pari_mutex();
1017     if (!locked)
1018       return 0;
1019 #ifdef HAVE_LIBPTHREAD
1020     delete pari_mutex_ptr;
1021     pari_mutex_ptr = 0;
1022 #endif
1023     return 1;
1024   }
1025   static const char _pari_unlock_s []="pari_unlock";
1026   static define_unary_function_eval (__pari_unlock,&_pari_unlock,_pari_unlock_s);
1027   define_unary_function_ptr5( at_pari_unlock ,alias_at_pari_unlock,&__pari_unlock,_QUOTE_ARGUMENTS,true);
1028 
cutstring(const std::string & s,int ncol)1029   static std::string cutstring(const std::string & s,int ncol){
1030     string res;
1031     if (ncol<20)
1032       ncol=20;
1033     int left=s.size(),pos=0,j;
1034     for (;left>ncol;pos+=j,left-=j){
1035       for (j=ncol;j>ncol/2;--j){
1036 	if (s[pos+j]==' ')
1037 	  break;
1038       }
1039       res=res+s.substr(pos,j)+'\n';
1040     }
1041     return res+s.substr(pos,pos+left);
1042   }
1043 
1044   // Help for g calling PARI online help
pari_help(const gen & g)1045   std::string pari_help(const gen & g){
1046     if (is_zero(g))
1047       return "Run pari() to export PARI functions.\n?pari(1) to ?pari(11) lists PARI functions by section\n?pari_functionname shows a short help on a function\nInside Xcas, Help->Manual->PARI-GP shows HTML help";
1048     string res;
1049     if (g.type==_INT_){
1050       int section=g.val;
1051       entree * ptr=functions_basic;
1052       for (;ptr->name;++ptr){
1053 	if (ptr->menu==section){
1054 	  res += ptr->name;
1055 	  res += " ";
1056 	}
1057       }
1058       return cutstring(res,70);
1059     }
1060     string gs;
1061     if (g.type==_FUNC)
1062       gs=g._FUNCptr->ptr()->s;
1063     else
1064       gs=gen2string(g);
1065     if (gs.size()>5 && gs.substr(0,5)=="pari_")
1066       gs=gs.substr(5,gs.size()-5);
1067     entree * ptr=functions_basic;
1068     for (;ptr->name;++ptr){
1069       if (ptr->name==gs){
1070 	res = ptr->help;
1071 	return cutstring(res,70);
1072       }
1073     }
1074     return "PARI function not found\nHelp syntax: ?pari(1),...,?pari(12) or ?pari_functionname";
1075   }
1076 
pari_polroots(const vecteur & p,vecteur & res,long prec,GIAC_CONTEXT)1077   bool pari_polroots(const vecteur & p,vecteur & res,long prec,GIAC_CONTEXT){
1078     if (check_pari_mutex())
1079       return false;
1080     gen tmp;
1081 #ifdef HAVE_LIBPTHREAD
1082     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
1083 #endif
1084     long av=get_pari_avma();
1085     GEN G=gen2GEN(change_subtype(p,_POLY1__VECT),vecteur(0),contextptr);
1086     if (debug_infolevel)
1087       CERR << "pari_polroots " << GEN2gen(G,vecteur(1,vx_var)) << '\n';
1088     G=roots(G,prec);
1089     tmp=GEN2gen(G,vecteur(0));
1090     avma=av;
1091 #ifdef HAVE_LIBPTHREAD
1092     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
1093     pthread_cleanup_pop(0);
1094 #endif
1095     if (tmp.type!=_VECT)
1096       return false;
1097     res=*tmp._VECTptr;
1098     return true;
1099   }
1100 
pari_polresultant(const gen & p,const gen & q,const vecteur & lv,gen & res,GIAC_CONTEXT)1101   bool pari_polresultant(const gen & p,const gen & q,const vecteur & lv,gen & res,GIAC_CONTEXT){
1102     if (check_pari_mutex())
1103       return false;
1104     gen tmp;
1105 #ifdef HAVE_LIBPTHREAD
1106     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
1107 #endif
1108     long av=get_pari_avma();
1109     void * save_pari_stack_limit = PARI_stack_limit;
1110     PARI_stack_limit=0;
1111     GEN P=gen2GEN(p,lv,contextptr);
1112     GEN Q=gen2GEN(q,lv,contextptr);
1113     GEN PQ=polresultant0(P,Q,-1,2);
1114     tmp=GEN2gen(PQ,lv);
1115     avma=av;
1116     PARI_stack_limit=save_pari_stack_limit;
1117 #ifdef HAVE_LIBPTHREAD
1118     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
1119     pthread_cleanup_pop(0);
1120 #endif
1121     res=tmp;
1122     return true;
1123   }
1124 
pari_nffactor(const gen & p,const gen & pmin,const vecteur & lv,gen & res,GIAC_CONTEXT)1125   bool pari_nffactor(const gen & p,const gen & pmin,const vecteur & lv,gen & res,GIAC_CONTEXT){
1126     if (check_pari_mutex())
1127       return false;
1128     gen tmp;
1129 #ifdef HAVE_LIBPTHREAD
1130     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
1131 #endif
1132     long av=get_pari_avma();
1133     void * save_pari_stack_limit = PARI_stack_limit;
1134     PARI_stack_limit=0;
1135     GEN P=gen2GEN(p,lv,contextptr);
1136     GEN Pmin=gen2GEN(pmin,lv,contextptr);
1137     int prec=decimal_digits(contextptr);
1138     if (prec<30)
1139       prec=30;
1140 #if 0
1141     GEN nf=nfinit0(Pmin,0,prec);
1142     tmp=GEN2gen(nffactor(nf,P),lv);
1143 #else
1144     tmp=GEN2gen(nffactor(Pmin,P),lv);
1145 #endif
1146     avma=av;
1147     PARI_stack_limit=save_pari_stack_limit;
1148 #ifdef HAVE_LIBPTHREAD
1149     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
1150     pthread_cleanup_pop(0);
1151 #endif
1152     res=tmp;
1153     return true;
1154   }
1155 
pari_galoisconj(const gen & g,vecteur & w,GIAC_CONTEXT)1156   bool pari_galoisconj(const gen & g,vecteur & w,GIAC_CONTEXT){
1157     if (check_pari_mutex())
1158       return false;
1159     gen res;
1160 #ifdef HAVE_LIBPTHREAD
1161     pthread_cleanup_push(pari_cleanup, (void *) pari_mutex_ptr);
1162 #endif
1163     res=in_pari(makesequence(string2gen("nfgaloisconj",false),g),contextptr);
1164 #ifdef HAVE_LIBPTHREAD
1165     if (pari_mutex_ptr) pthread_mutex_unlock(pari_mutex_ptr);
1166     pthread_cleanup_pop(0);
1167 #endif
1168     if (res.type!=_VECT)
1169       return false;
1170     w=*res._VECTptr;
1171     gen gp=_symb2poly(makesequence(g,vx_var),contextptr);
1172     for (int i=0;i<int(w.size());++i){
1173       gen tmp=w[i];
1174       tmp=_symb2poly(makesequence(tmp,vx_var),contextptr);
1175       gen d=1;
1176       if (tmp.type==_VECT)
1177 	lcmdeno(*tmp._VECTptr,d,contextptr);
1178       w[i]=symb_rootof(tmp,gp,contextptr)/d;
1179     }
1180     return true;
1181   }
1182 
1183 #ifndef NO_NAMESPACE_GIAC
1184 }
1185 #endif // ndef NO_NAMESPACE_GIAC
1186 
1187 #else
1188 // ! HAVE_LIBPARI
1189 #ifndef NO_NAMESPACE_GIAC
1190 namespace giac {
1191 #endif // ndef NO_NAMESPACE_GIAC
1192 
pari_error()1193   static gen pari_error(){
1194     return undeferr(gettext("Not implemented, please recompile giac with PARI"));
1195   }
1196 
pari_isprime(const gen & e,int certif)1197   gen pari_isprime(const gen & e,int certif){
1198     return string2gen("please recompile giac with PARI",false);
1199   }
pari_ifactor(const gen & e)1200   std::string pari_ifactor(const gen & e){
1201     return "please recompile giac with PARI";
1202   }
1203 
pari_ffinit(const gen & p,int n)1204   gen pari_ffinit(const gen & p,int n){
1205     return string2gen("please recompile giac with PARI",false);
1206   }
1207 
pari_gamma(const gen & e)1208   gen pari_gamma(const gen & e){
1209     return pari_error();
1210   }
pari_zeta(const gen & e)1211   gen pari_zeta(const gen & e){
1212     return pari_error();
1213   }
pari_psi(const gen & e)1214   gen pari_psi(const gen & e){
1215     return pari_error();
1216   }
1217 
pari_lift_combine(const vecteur & a,const std::vector<vecteur> & factmod,gen & modulo,std::vector<vecteur> & res)1218   bool pari_lift_combine(const vecteur& a, const std::vector<vecteur>& factmod,
1219 			 gen& modulo, std::vector<vecteur>& res){
1220     vecteur tmp(1,pari_error());
1221     res=std::vector<vecteur>(1,tmp);
1222     return false;
1223   }
1224 
pari_galoisconj(const gen & g,vecteur & w,GIAC_CONTEXT)1225   bool pari_galoisconj(const gen & g,vecteur & w,GIAC_CONTEXT){
1226     return false;
1227   }
1228 
_pari(const gen & args,GIAC_CONTEXT)1229   gen _pari(const gen & args,GIAC_CONTEXT){
1230     if ( args.type==_STRNG && args.subtype==-1) return  args;
1231     return pari_error();
1232   }
1233 
pari_help(const gen & g)1234   std::string pari_help(const gen & g){
1235     return "please recompile giac with PARI";
1236   }
1237 
GEN2gen(const GEN & G,const vecteur & vars)1238   gen GEN2gen(const GEN & G,const vecteur & vars){
1239     return gensizeerr("please recompile giac with PARI");
1240   }
gen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT)1241   GEN gen2GEN(const gen & e,const vecteur & vars,GIAC_CONTEXT){
1242     return 0;
1243   }
1244 
pari_polroots(const vecteur & p,vecteur & res,long prec,GIAC_CONTEXT)1245   bool pari_polroots(const vecteur & p,vecteur & res,long prec,GIAC_CONTEXT){
1246     return false;
1247   }
pari_polresultant(const gen & p,const gen & q,const vecteur & lv,gen & res,GIAC_CONTEXT)1248   bool pari_polresultant(const gen & p,const gen & q,const vecteur & lv,gen & res,GIAC_CONTEXT){
1249     return false;
1250   }
pari_nffactor(const gen & p,const gen & pmin,const vecteur & lv,gen & res,GIAC_CONTEXT)1251   bool pari_nffactor(const gen & p,const gen & pmin,const vecteur & lv,gen & res,GIAC_CONTEXT){
1252     return false;
1253   }
1254   static const char _pari_s []="pari";
1255   static define_unary_function_eval (__pari,&_pari,_pari_s);
1256   define_unary_function_ptr5( at_pari ,alias_at_pari,&__pari,_QUOTE_ARGUMENTS,true);
1257 
1258   static const char _pari_unlock_s []="pari_unlock";
1259   static define_unary_function_eval (__pari_unlock,&_pari,_pari_unlock_s);
1260   define_unary_function_ptr5( at_pari_unlock ,alias_at_pari_unlock,&__pari_unlock,_QUOTE_ARGUMENTS,true);
1261 
1262 
1263 #ifndef NO_NAMESPACE_GIAC
1264 }
1265 #endif // ndef NO_NAMESPACE_GIAC
1266 
1267 #endif // HAVE_LIBPARI
1268