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