=================================================================== RCS file: /home/cvs/OpenXM/src/ox_pari/ox_pari.c,v retrieving revision 1.8 retrieving revision 1.16 diff -u -p -r1.8 -r1.16 --- OpenXM/src/ox_pari/ox_pari.c 2015/08/20 01:38:34 1.8 +++ OpenXM/src/ox_pari/ox_pari.c 2018/03/27 07:05:17 1.16 @@ -1,49 +1,15 @@ -/* $OpenXM: OpenXM/src/ox_pari/ox_pari.c,v 1.7 2015/08/18 05:04:35 noro Exp $ */ +/* $OpenXM: OpenXM/src/ox_pari/ox_pari.c,v 1.15 2016/09/23 07:03:29 noro Exp $ */ -#include -#include -#include -#include "pari/pari.h" -#include "pari/paripriv.h" -#include "gmp.h" -#include "gmp-impl.h" -#include "mpfr.h" -#include "ox_toolkit.h" +#include "ox_pari.h" + OXFILE *fd_rw; -#define MPFR_PREC(x) ((x)->_mpfr_prec) -#define MPFR_EXP(x) ((x)->_mpfr_exp) -#define MPFR_MANT(x) ((x)->_mpfr_d) -#define MPFR_LAST_LIMB(x) ((MPFR_PREC (x) - 1) / GMP_NUMB_BITS) -#define MPFR_LIMB_SIZE(x) (MPFR_LAST_LIMB (x) + 1) - static int stack_size = 0; static int stack_pointer = 0; static cmo **stack = NULL; extern int debug_print; long paristack=10000000; -void init_pari(void); -cmo *GEN_to_cmo(GEN z); -cmo_zz *GEN_to_cmo_zz(GEN z); -cmo_qq *GEN_to_cmo_qq(GEN z); -cmo_bf *GEN_to_cmo_bf(GEN z); -cmo_list *GEN_to_cmo_list(GEN z); -cmo_complex *GEN_to_cmo_cmo_complex(GEN z); -cmo_polynomial_in_one_variable *GEN_to_cmo_up(GEN z); -cmo_recursive_polynomial *GEN_to_cmo_rp(GEN z); - -GEN cmo_to_GEN(cmo *c); -GEN cmo_int32_to_GEN(cmo_int32 *c); -GEN cmo_zz_to_GEN(cmo_zz *c); -GEN cmo_qq_to_GEN(cmo_qq *c); -GEN cmo_bf_to_GEN(cmo_bf *c); -GEN cmo_list_to_GEN(cmo_list *c); -GEN cmo_rp_to_GEN(cmo_recursive_polynomial *c); -GEN cmo_up_to_GEN(cmo_polynomial_in_one_variable *c); -GEN cmo_complex_to_GEN(cmo_complex *c); - - #define INIT_S_SIZE 2048 #define EXT_S_SIZE 2048 @@ -119,7 +85,12 @@ void pops(int n) int sm_mathcap() { - mathcap_init(OX_PARI_VERSION, ID_STRING, "ox_pari", NULL, NULL); +#if 0 + char *opts[] = {"no_ox_reset", NULL}; + mathcap_init2(OX_PARI_VERSION, ID_STRING, "ox_pari", NULL, NULL, opts); +#else + mathcap_init2(OX_PARI_VERSION, ID_STRING, "ox_pari", NULL, NULL, NULL); +#endif push((cmo*)oxf_cmo_mathcap(fd_rw)); return 0; } @@ -137,7 +108,7 @@ int sm_popCMO() cmo_error2 *make_error2(char *message) { - return (cmo_error2 *) new_cmo_string(message); + return new_cmo_error2((cmo *)new_cmo_string(message)); } int get_i() @@ -162,392 +133,6 @@ char *get_str() return ""; } -GEN cmo_int32_to_GEN(cmo_int32 *c) -{ - GEN z; - int i,sgn; - - i = c->i; - if ( !i ) return gen_0; - z = cgeti(3); - sgn = 1; - if ( i < 0 ) { - i = -i; - sgn = -1; - } - z[2] = i; - setsigne(z,sgn); - setlgefint(z,lg(z)); - return z; -} - -GEN cmo_zz_to_GEN(cmo_zz *c) -{ - mpz_ptr mpz; - GEN z; - long *ptr; - int j,sgn,len; - - mpz = c->mpz; - sgn = mpz_sgn(mpz); - len = ABSIZ(mpz); - ptr = (long *)PTR(mpz); - z = cgeti(len+2); - for ( j = 0; j < len; j++ ) - z[len-j+1] = ptr[j]; - setsigne(z,sgn); - setlgefint(z,lg(z)); - return z; -} - -GEN cmo_qq_to_GEN(cmo_qq *c) -{ - GEN z,nm,den; - - z = cgetg(3,t_FRAC); - nm = cmo_zz_to_GEN(new_cmo_zz_set_mpz(mpq_numref(c->mpq))); - den = cmo_zz_to_GEN(new_cmo_zz_set_mpz(mpq_denref(c->mpq))); - z[1] = (long)nm; - z[2] = (long)den; - return z; -} - -GEN cmo_bf_to_GEN(cmo_bf *c) -{ - mpfr_ptr mpfr; - GEN z; - int sgn,len,j; - long exp; - long *ptr; - - mpfr = c->mpfr; - sgn = MPFR_SIGN(mpfr); - exp = MPFR_EXP(mpfr)-1; - len = MPFR_LIMB_SIZE(mpfr); - ptr = (long *)MPFR_MANT(mpfr); - z = cgetr(len+2); - for ( j = 0; j < len; j++ ) - z[len-j+1] = ptr[j]; - z[1] = evalsigne(sgn)|evalexpo(exp); - setsigne(z,sgn); - return z; -} - -/* list->vector */ - -GEN cmo_list_to_GEN(cmo_list *c) -{ - GEN z; - int i; - cell *cell; - - z = cgetg(c->length+1,t_VEC); - for ( i = 0, cell = c->head->next; cell != c->head; cell = cell->next, i++ ) { - z[i+1] = (long)cmo_to_GEN(cell->cmo); - } - return z; -} - -GEN cmo_complex_to_GEN(cmo_complex *c) -{ - GEN z; - - z = cgetg(3,t_COMPLEX); - z[1] = (long)cmo_to_GEN(c->re); - z[2] = (long)cmo_to_GEN(c->im); - return z; -} - -GEN cmo_up_to_GEN(cmo_polynomial_in_one_variable *c) -{ - GEN z; - int d,i; - cell *cell; - - d = c->head->next->exp; - z = cgetg(d+3,t_POL); - setsigne(z,1); - setvarn(z,c->var); - setlgef(z,d+3); - for ( i = 2; i <= d+2; i++ ) - z[i] = (long)gen_0; - for ( cell = c->head->next; cell != c->head; cell = cell->next ) { - z[2+cell->exp] = (long)cmo_to_GEN(cell->cmo); - } - return z; -} - -cmo_list *current_ringdef; - -void register_variables(cmo_list *ringdef) -{ - current_ringdef = ringdef; -} - -GEN cmo_rp_to_GEN(cmo_recursive_polynomial *c) -{ - register_variables(c->ringdef); - switch ( c->coef->tag ) { - case CMO_ZERO: - case CMO_NULL: - return gen_0; - case CMO_INT32: - return cmo_int32_to_GEN((cmo_int32 *)c->coef); - case CMO_ZZ: - return cmo_zz_to_GEN((cmo_zz *)c->coef); - case CMO_QQ: - return cmo_qq_to_GEN((cmo_qq *)c->coef); - case CMO_POLYNOMIAL_IN_ONE_VARIABLE: - return cmo_up_to_GEN((cmo_polynomial_in_one_variable *)c->coef); - default: - return 0; - } -} - -cmo_zz *GEN_to_cmo_zz(GEN z) -{ - cmo_zz *c; - - c = new_cmo_zz(); - mpz_import(c->mpz,lgef(z)-2,1,sizeof(long),0,0,&z[2]); - if ( signe(z) < 0 ) - mpz_neg(c->mpz,c->mpz); - return c; -} - -cmo_qq *GEN_to_cmo_qq(GEN z) -{ - cmo_qq *c; - GEN num,den; - - num = (GEN)z[1]; - den = (GEN)z[2]; - c = new_cmo_qq(); - mpz_import(mpq_numref(c->mpq),lgef(num)-2,1,sizeof(long),0,0,&num[2]); - mpz_import(mpq_denref(c->mpq),lgef(num)-2,1,sizeof(long),0,0,&den[2]); - if ( signe(num)*signe(den) < 0 ) - mpz_neg(mpq_numref(c->mpq),mpq_numref(c->mpq)); - return c; -} - - -cmo_bf *GEN_to_cmo_bf(GEN z) -{ - cmo_bf *c; - int len,prec,j; - long *ptr; - - c = new_cmo_bf(); - len = lg(z)-2; - prec = len*sizeof(long)*8; - mpfr_init2(c->mpfr,prec); - ptr = (long *)MPFR_MANT(c->mpfr); - for ( j = 0; j < len; j++ ) - ptr[j] = z[len-j+1]; - MPFR_EXP(c->mpfr) = (long long)(expo(z)+1); - MPFR_SIGN(c->mpfr) = gsigne(z); - return c; -} - - -cmo_list *GEN_to_cmo_list(GEN z) -{ - cmo_list *c; - cmo *ob; - int i,len; - - c = new_cmo_list(); - len = lg(z)-1; - for ( i = 1; i <= len; i++ ) { - ob = GEN_to_cmo((GEN)z[i]); - c = list_append(c,ob); - } - return c; -} - -cmo_complex *GEN_to_cmo_complex(GEN z) -{ - cmo_complex *c; - - c = new_cmo_complex(); - c->re = GEN_to_cmo((GEN)z[1]); - c->im = GEN_to_cmo((GEN)z[2]); - return c; -} - -cmo_polynomial_in_one_variable *GEN_to_cmo_up(GEN z) -{ - cmo_polynomial_in_one_variable *c; - int i; - cmo *coef; - - c = new_cmo_polynomial_in_one_variable(varn(z)); - for ( i = lg(z)-1; i >= 2; i-- ) - if ( (GEN)z[i] != gen_0 ) { - coef = GEN_to_cmo((GEN)z[i]); - list_append_monomial((cmo_list *)c, coef, i-2); - } - return c; -} - -cmo_recursive_polynomial *GEN_to_cmo_rp(GEN z) -{ - cmo_recursive_polynomial *c; - - if ( !signe(z) ) return (cmo_recursive_polynomial *)new_cmo_zero(); - c = new_cmo_recursive_polynomial(current_ringdef,(cmo *)GEN_to_cmo_up(z)); - return c; -} - -GEN cmo_to_GEN(cmo *c) -{ - switch ( c->tag ) { - case CMO_ZERO: - case CMO_NULL: - return gen_0; - case CMO_ZZ: /* int */ - return cmo_zz_to_GEN((cmo_zz *)c); - case CMO_IEEE_DOUBLE_FLOAT: - return dbltor(((cmo_double *)c)->d); - case CMO_BIGFLOAT: /* bigfloat */ - return cmo_bf_to_GEN((cmo_bf *)c); - case CMO_LIST: - return cmo_list_to_GEN((cmo_list *)c); - case CMO_RECURSIVE_POLYNOMIAL: - return cmo_rp_to_GEN((cmo_recursive_polynomial *)c); - case CMO_POLYNOMIAL_IN_ONE_VARIABLE: - return cmo_up_to_GEN((cmo_polynomial_in_one_variable *)c); - default: - return 0; - } -} - -cmo *GEN_to_cmo(GEN z) -{ - char buf[BUFSIZ]; - - if ( gcmp0(z) ) - return new_cmo_zero(); - switch ( typ(z) ) { - case t_INT: /* int */ - return (cmo *)GEN_to_cmo_zz(z); - case t_REAL: /* bigfloat */ - return (cmo *)GEN_to_cmo_bf(z); - case t_FRAC: /* rational number */ - return (cmo *)GEN_to_cmo_qq(z); - case t_COMPLEX: /* complex */ - return (cmo *)GEN_to_cmo_complex(z); - case t_POL: - return (cmo *)GEN_to_cmo_rp(z); - case t_VEC: case t_COL: /* vector */ - return (cmo *)GEN_to_cmo_list(z); - case t_MAT: /* matrix */ - return (cmo *)GEN_to_cmo_list(shallowtrans(z)); - default: - sprintf(buf,"GEN_to_cmo : unsupported type=%d",(int)typ(z)); - return (cmo *)make_error2(buf); - } -} -/* type=1 : num/poly arg, type=2 : matrix arg */ - -struct parif { - char *name; - GEN (*f)(); - int type; -} parif_tab[] = { -/* (ulong)allocatemoremem(ulong) */ - {"allocatemem",(GEN (*)())allocatemoremem,0}, -/* num/num */ - {"abs",gabs,1}, - {"erfc",gerfc,1}, - {"arg",garg,1}, - {"isqrt",racine,1}, - {"gamma",ggamma,1}, - {"zeta",gzeta,1}, - {"floor",gfloor,1}, - {"frac",gfrac,1}, - {"imag",gimag,1}, - {"conj",gconj,1}, - {"ceil",gceil,1}, - {"isprime",gisprime,2}, - {"bigomega",gbigomega,1}, - {"denom",denom,1}, - {"numer",numer,1}, - {"lngamma",glngamma,1}, - {"logagm",glogagm,1}, - {"classno",classno,1}, - {"classno2",classno2,1}, - {"dilog",dilog,1}, - {"disc",discsr,1}, - {"discf",discf,1}, - {"nextprime",nextprime,1}, - {"eintg1",eint1,1}, - {"eta",eta,1}, - {"issqfree",gissquarefree,1}, - {"issquare",gcarreparfait,1}, - {"gamh",ggamd,1}, - {"hclassno",classno3,1}, - - /* num/array */ - {"binary",binaire,1}, - {"factorint",factorint,2}, - {"factor",Z_factor,1}, - {"cf",gcf,1}, - {"divisors",divisors,1}, - {"smallfact",smallfact,1}, - - /* poly/poly */ - {"centerlift",centerlift,1}, - {"content",content,1}, - - /* poly/array */ - {"galois",galois,1}, - {"roots",roots,1}, - - /* mat/mat */ - {"adj",adj,2}, - {"lll",lll,2}, - {"lllgen",lllgen,2}, - {"lllgram",lllgram,2}, - {"lllgramgen",lllgramgen,2}, - {"lllgramint",lllgramint,2}, - {"lllgramkerim",lllgramkerim,2}, - {"lllgramkerimgen",lllgramkerimgen,2}, - {"lllint",lllint,2}, - {"lllkerim",lllkerim,2}, - {"lllkerimgen",lllkerimgen,2}, - {"trans",gtrans,2}, - {"eigen",eigen,2}, - {"hermite",hnf,2}, - {"mat",gtomat,2}, - {"matrixqz2",matrixqz2,2}, - {"matrixqz3",matrixqz3,2}, - {"hess",hess,2}, - {"ker",ker,2}, - {"keri",keri,2}, - {"kerint",kerint,2}, - {"kerintg1",kerint1,2}, - - /* mat/poly */ - {"det",det,2}, - {"det2",det2,2}, - -}; - -#define PARI_MAX_AC 64 - -struct parif *search_parif(char *name) -{ - int tablen,i; - - tablen = sizeof(parif_tab)/sizeof(struct parif); - for ( i = 0; i < tablen; i++ ) { - if ( !strcmp(parif_tab[i].name,name) ) - return &parif_tab[i]; - } - return 0; -} - int ismatrix(GEN z) { int len,col,i; @@ -571,9 +156,18 @@ int sm_executeFunction() struct parif *parif; unsigned long prec; char buf[BUFSIZ]; + int status; + char *err; - if ( setjmp(GP_DATA->env) ) { - sprintf(buf,"sm_executeFunction : an error occured in PARI."); + if ( (status = setjmp(GP_DATA->env)) != 0 ) { + err = errmessage[status]; + if ( status == errpile ) { + sprintf(buf,"%s\nIncrease PARI stack by pari(allocatemem,size).",err); + init_pari(); + } else if ( strlen(err) != 0 ) + sprintf(buf,"An error occured in PARI :%s",err); + else + sprintf(buf,"An error occured in PARI."); push((cmo*)make_error2(buf)); return -1; } @@ -611,7 +205,7 @@ int sm_executeFunction() ret = (cmo *)new_cmo_int32(a); push(ret); return 0; - } else if ( parif->type == 1 || parif->type == 2 ) { + } else if ( parif->type == 1 ) { /* one number/poly/matrix argument possibly with prec */ av0 = avma; z = cmo_to_GEN(av[0]); @@ -630,6 +224,24 @@ int sm_executeFunction() avma = av0; push(ret); return 0; + } else if ( parif->type == 2 ) { + /* one number/poly/matrix argument with flag=0 */ + av0 = avma; + z = cmo_to_GEN(av[0]); + if ( ismatrix(z) ) { + int i,len; + len = lg(z); + for ( i = 1; i < len; i++ ) + settyp(z[i],t_COL); + settyp(z,t_MAT); + z = shallowtrans(z); + } + printf("input : "); output(z); + m = (*parif->f)(z,0); + ret = GEN_to_cmo(m); + avma = av0; + push(ret); + return 0; } else { sprintf(buf,"%s : not implemented",func->s); push((cmo*)make_error2(buf)); @@ -653,6 +265,9 @@ int receive_and_execute_sm_command() case SM_setMathCap: pop(); break; + case SM_shutdown: + exit(0); + break; default: printf("receive_and_execute_sm_command : code=%d\n",code);fflush(stdout); break; @@ -680,17 +295,33 @@ int receive() return 0; } +jmp_buf ox_env; + +void usr1_handler(int sig) +{ + longjmp(ox_env,1); +} + int main() { - init_gc(); - ox_stderr_init(stderr); - initialize_stack(); - init_pari(); + if ( setjmp(ox_env) ) { + fprintf(stderr,"resetting libpari and sending OX_SYNC_BALL..."); + initialize_stack(); + init_pari(); + send_ox_tag(fd_rw,OX_SYNC_BALL); + fprintf(stderr,"done\n"); + } else { + init_gc(); + ox_stderr_init(stderr); + initialize_stack(); + init_pari(); + + fprintf(stderr,"ox_pari\n"); - fprintf(stderr,"ox_pari\n"); - - fd_rw = oxf_open(3); - oxf_determine_byteorder_server(fd_rw); + fd_rw = oxf_open(3); + oxf_determine_byteorder_server(fd_rw); + } + signal(SIGUSR1,usr1_handler); while(1){ receive();