=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/parif.c,v retrieving revision 1.9 retrieving revision 1.17 diff -u -p -r1.9 -r1.17 --- OpenXM_contrib2/asir2000/builtin/parif.c 2001/10/03 01:47:30 1.9 +++ OpenXM_contrib2/asir2000/builtin/parif.c 2011/08/10 04:51:57 1.17 @@ -45,21 +45,27 @@ * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. * - * $OpenXM: OpenXM_contrib2/asir2000/builtin/parif.c,v 1.8 2001/07/05 09:26:43 noro Exp $ + * $OpenXM: OpenXM_contrib2/asir2000/builtin/parif.c,v 1.16 2009/03/27 14:42:29 ohara Exp $ */ #include "ca.h" #include "parse.h" -#if PARI +#if defined(PARI) #include "genpari.h" +// PARI_VERSION(2,2,8) == 131592 +#if PARI_VERSION_CODE >= 131592 +#define mppgcd(a,b) (gcdii((a),(b))) +#endif -extern long prec; +long get_pariprec(); +void set_pariprec(long p); void patori(GEN,Obj *); void patori_i(GEN,N *); void ritopa(Obj,GEN *); void ritopa_i(N,int,GEN *); +void Ptodouble(); void Peval(),Psetprec(),p_pi(),p_e(),p_mul(),p_gcd(); void asir_cgiv(GEN); @@ -70,13 +76,14 @@ void Psetprecword(); struct ftab pari_tab[] = { {"eval",Peval,-2}, {"setprec",Psetprec,-1}, + {"todouble",Ptodouble,1}, #if defined(INTERVAL) || 1 {"setprecword",Psetprecword,-1}, #endif {0,0,0}, }; -#define MKPREC(a,i,b) (argc(a)==(i)?mkprec(QTOS((Q)(b))):prec) +#define MKPREC(a,i,b) (argc(a)==(i)?mkprec(QTOS((Q)(b))):get_pariprec()) #define CALLPARI1(f,a,p,r)\ ritopa((Obj)a,&_pt1_); _pt2_ = f(_pt1_,p); patori(_pt2_,r); asir_cgiv(_pt2_); asir_cgiv(_pt1_) @@ -93,15 +100,12 @@ void f(ar,rp) NODE ar; Obj *rp;\ { GEN _pt1_,_pt2_,_pt3_; CALLPARI2(pf,ARG0(ar),ARG1(ar),MKPREC(ar,3,ARG2(ar)),rp); } #if defined(INTERVAL) -#define PREC_CONV pariK1 -#else -#if defined(LONG_IS_32BIT) +#define PREC_CONV (0.103810253/(BYTES_IN_LONG/4)) /* log(10)/(SL*log(2)) */ +#elif SIZEOF_LONG == 4 #define PREC_CONV 0.103810253 -#endif -#if defined(LONG_IS_64BIT) +#elif SIZEOF_LONG == 8 #define PREC_CONV 0.051905126 #endif -#endif /* XXX : we should be more careful when we free PARI pointers. */ @@ -117,10 +121,44 @@ GEN ptr; mkprec(p) int p; { - if ( p > 0 ) - return (int)(p*PREC_CONV+3); + if ( p <= 0 ) + p = 1; + return (int)(p*PREC_CONV+3); } +void Ptodouble(arg,rp) +NODE arg; +Num *rp; +{ + double r,i; + Real real,imag; + Num num; + + asir_assert(ARG0(arg),O_N,"todouble"); + num = (Num)ARG0(arg); + if ( !num ) { + *rp = 0; + return; + } + switch ( NID(num) ) { + case N_R: case N_Q: case N_B: + r = ToReal(num); + MKReal(r,real); + *rp = (Num)real; + break; + case N_C: + r = ToReal(((C)num)->r); + i = ToReal(((C)num)->i); + MKReal(r,real); + MKReal(i,imag); + reimtocplx((Num)real,(Num)imag,rp); + break; + default: + *rp = num; + break; + } +} + void Peval(arg,rp) NODE arg; Obj *rp; @@ -135,6 +173,7 @@ Obj *rp; { int p; Q q; + long prec = get_pariprec(); #if defined(INTERVAL) || 1 p = (int)((prec-2)/PREC_CONV); STOQ(p,q); *rp = (Obj)q; @@ -151,6 +190,7 @@ Obj *rp; prec = mkprec(QTOS((Q)ARG0(arg))); } #endif + set_pariprec(prec); } #if defined(INTERVAL) || 1 @@ -160,6 +200,7 @@ Obj *rp; { int p; Q q; + long prec = get_pariprec(); p = (int)((prec-2)); STOQ(p,q); *rp = (Obj)q; if ( arg ) { @@ -169,6 +210,7 @@ Obj *rp; prec = p + 2; } } + set_pariprec(prec); } #endif @@ -236,6 +278,8 @@ NODE arg; if ( !f->f.binf ) { sprintf(buf,"pari : %s undefined.",f->name); error(buf); + /* NOTREACHED */ + return 0; } switch ( f->type ) { case 0: /* in/out : integer */ @@ -243,10 +287,12 @@ NODE arg; if ( ac > 2 ) { fprintf(stderr,"argument mismatch in %s()\n",NAME(f)); error(""); + /* NOTREACHED */ + return 0; } intarg = !ac ? 0 : QTOS((Q)ARG0(arg)); dmy = (GEN (*)())f->f.binf; - ret = (*dmy)(intarg); + ret = (int)(*dmy)(intarg); STOQ(ret,q); return (pointer)q; @@ -255,6 +301,8 @@ NODE arg; if ( !ac || ( ac > 2 ) ) { fprintf(stderr,"argument mismatch in %s()\n",NAME(f)); error(""); + /* NOTREACHED */ + return 0; } ltop = avma; ritopa((Obj)ARG0(arg),&a); @@ -269,6 +317,8 @@ NODE arg; if ( !ac || ( ac > 2 ) ) { fprintf(stderr,"argument mismatch in %s()\n",NAME(f)); error(""); + /* NOTREACHED */ + return 0; } if ( ac == 1 ) opt = 0; @@ -284,6 +334,8 @@ NODE arg; default: error("evalparif : not implemented yet."); + /* NOTREACHED */ + return 0; } } @@ -300,7 +352,7 @@ struct pariftab { */ struct pariftab pariftab[] = { -{"allocatemem",allocatemoremem,0}, +{"allocatemem",(GEN(*)())allocatemoremem,0}, {"abs",(GEN (*)())gabs,1}, {"adj",adj,1}, {"arg",garg,1}, @@ -340,7 +392,7 @@ struct pariftab pariftab[] = { {"indsort",indexsort,1}, {"initalg",initalg,1}, {"isfund",gisfundamental,1}, -{"isprime",gisprime,1}, +{"isprime",gisprime,2}, {"ispsp",gispsp,1}, {"isqrt",racine,1}, {"issqfree",gissquarefree,1}, @@ -351,16 +403,13 @@ struct pariftab pariftab[] = { {"keri",keri,1}, {"kerint",kerint,1}, {"kerintg1",kerint1,1}, -{"kerint2",kerint2,1}, {"length",(GEN(*)())glength,1}, {"lexsort",lexsort,1}, {"lift",lift,1}, {"lindep",lindep,1}, {"lll",lll,1}, -{"lllg1",lll1,1}, {"lllgen",lllgen,1}, {"lllgram",lllgram,1}, -{"lllgramg1",lllgram1,1}, {"lllgramgen",lllgramgen,1}, {"lllgramint",lllgramint,1}, {"lllgramkerim",lllgramkerim,1}, @@ -368,7 +417,6 @@ struct pariftab pariftab[] = { {"lllint",lllint,1}, {"lllkerim",lllkerim,1}, {"lllkerimgen",lllkerimgen,1}, -{"lllrat",lllrat,1}, {"lngamma",glngamma,1}, {"logagm",glogagm,1}, {"mat",gtomat,1},