=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/bfaux.c,v retrieving revision 1.1 retrieving revision 1.4 diff -u -p -r1.1 -r1.4 --- OpenXM_contrib2/asir2000/builtin/bfaux.c 2015/08/04 06:55:02 1.1 +++ OpenXM_contrib2/asir2000/builtin/bfaux.c 2015/08/06 23:41:52 1.4 @@ -1,16 +1,20 @@ -/* $OpenXM$ */ +/* $OpenXM: OpenXM_contrib2/asir2000/builtin/bfaux.c,v 1.3 2015/08/05 01:23:02 noro Exp $ */ #include "ca.h" #include "parse.h" -void Peval(), Psetprec(), Ptodouble(); +void Peval(), Psetprec(), Psetbprec(), Ptodouble(), Psetround(); struct ftab bf_tab[] = { {"eval",Peval,-2}, {"setprec",Psetprec,-1}, + {"setbprec",Psetbprec,-1}, + {"setround",Psetround,-1}, {"todouble",Ptodouble,1}, {0,0,0}, }; +int mpfr_roundmode = MPFR_RNDN; + void Ptodouble(NODE arg,Num *rp) { double r,i; @@ -42,35 +46,56 @@ void Ptodouble(NODE arg,Num *rp) } } -void Peval(arg,rp) -NODE arg; -Obj *rp; +void Peval(NODE arg,Obj *rp) { int prec; asir_assert(ARG0(arg),O_R,"eval"); if ( argc(arg) == 2 ) { - prec = QTOS((Q)ARG1(arg)); + prec = QTOS((Q)ARG1(arg))*3.32193; if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN; else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX; } else prec = 0; - evalr(CO,(Obj)ARG0(arg),0,rp); + evalr(CO,(Obj)ARG0(arg),prec,rp); } -/* bit precision */ +/* set/get decimal precision */ void Psetprec(NODE arg,Obj *rp) { int p; Q q; + int prec,dprec; + + prec = mpfr_get_default_prec(); + /* decimal precision */ + dprec = prec*0.30103; + STOQ(dprec,q); *rp = (Obj)q; + if ( arg ) { + asir_assert(ARG0(arg),O_N,"setprec"); + prec = QTOS((Q)ARG0(arg))*3.32193; + if ( p > 0 ) + prec = p; + } + if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN; + else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX; + mpfr_set_default_prec(prec); +} - long prec = mpfr_get_default_prec(); +/* set/get bit precision */ - STOQ(prec,q); *rp = (Obj)q; +void Psetbprec(NODE arg,Obj *rp) +{ + int p; + Q q; + int prec; + + prec = mpfr_get_default_prec(); + STOQ(prec,q); *rp = (Obj)q; if ( arg ) { - asir_assert(ARG0(arg),O_N,"setprec"); - p = QTOS((Q)ARG0(arg)); + asir_assert(ARG0(arg),O_N,"setbprec"); + prec = QTOS((Q)ARG0(arg)); if ( p > 0 ) prec = p; } @@ -79,6 +104,43 @@ void Psetprec(NODE arg,Obj *rp) mpfr_set_default_prec(prec); } +void Psetround(NODE arg,Q *rp) +{ + int round; + + STOQ(mpfr_roundmode,*rp); + if ( arg ) { + asir_assert(ARG0(arg),O_N,"setround"); + round = QTOS((Q)ARG0(arg)); + switch ( round ) { + case 0: + mpfr_roundmode = MPFR_RNDN; + break; + case 1: + mpfr_roundmode = MPFR_RNDZ; + break; + case 2: + mpfr_roundmode = MPFR_RNDU; + break; + case 3: + mpfr_roundmode = MPFR_RNDD; + break; + case 4: + mpfr_roundmode = MPFR_RNDA; + break; + case 5: + mpfr_roundmode = MPFR_RNDF; + break; + case 6: + mpfr_roundmode = MPFR_RNDNA; + break; + default: + error("setround : invalid rounding mode"); + break; + } + } +} + Num tobf(Num a,int prec); void mp_pi(NODE arg,BF *rp) @@ -89,7 +151,7 @@ void mp_pi(NODE arg,BF *rp) prec = arg ? QTOS((Q)ARG0(arg)) : 0; NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_const_pi(r->body,MPFR_RNDN); + mpfr_const_pi(r->body,mpfr_roundmode); *rp = r; } @@ -103,8 +165,8 @@ void mp_e(NODE arg,BF *rp) NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); mpfr_init(one); - mpfr_set_ui(one,1,MPFR_RNDN); - mpfr_exp(r->body,one,MPFR_RNDN); + mpfr_set_ui(one,1,mpfr_roundmode); + mpfr_exp(r->body,one,mpfr_roundmode); *rp = r; } @@ -118,7 +180,7 @@ void mp_sin(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_sin(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_sin(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -132,7 +194,7 @@ void mp_cos(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_cos(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_cos(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -146,7 +208,7 @@ void mp_tan(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_tan(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_tan(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -160,7 +222,7 @@ void mp_asin(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_asin(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_asin(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } void mp_acos(NODE arg,BF *rp) @@ -173,7 +235,7 @@ void mp_acos(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_acos(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_acos(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } void mp_atan(NODE arg,BF *rp) @@ -186,7 +248,7 @@ void mp_atan(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_atan(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_atan(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -200,7 +262,7 @@ void mp_sinh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_sinh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_sinh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -214,7 +276,7 @@ void mp_cosh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_cosh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_cosh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -228,7 +290,7 @@ void mp_tanh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_tanh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_tanh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -242,7 +304,7 @@ void mp_asinh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_asinh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_asinh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } void mp_acosh(NODE arg,BF *rp) @@ -255,7 +317,7 @@ void mp_acosh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_acosh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_acosh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } void mp_atanh(NODE arg,BF *rp) @@ -268,7 +330,7 @@ void mp_atanh(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_atanh(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_atanh(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -282,7 +344,7 @@ void mp_exp(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_exp(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_exp(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -296,7 +358,7 @@ void mp_log(NODE arg,BF *rp) a = tobf(ARG0(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_log(r->body,((BF)a)->body,MPFR_RNDN); + mpfr_log(r->body,((BF)a)->body,mpfr_roundmode); *rp = r; } @@ -311,6 +373,6 @@ void mp_pow(NODE arg,BF *rp) e = tobf(ARG1(arg),prec); NEWBF(r); prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body); - mpfr_pow(r->body,((BF)a)->body,((BF)e)->body,MPFR_RNDN); + mpfr_pow(r->body,((BF)a)->body,((BF)e)->body,mpfr_roundmode); *rp = r; }