[BACK]Return to bfaux.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / builtin

Diff for /OpenXM_contrib2/asir2000/builtin/bfaux.c between version 1.1 and 1.8

version 1.1, 2015/08/04 06:55:02 version 1.8, 2015/08/16 03:12:09
Line 1 
Line 1 
 /* $OpenXM$ */  /* $OpenXM: OpenXM_contrib2/asir2000/builtin/bfaux.c,v 1.7 2015/08/07 08:00:30 takayama Exp $ */
 #include "ca.h"  #include "ca.h"
 #include "parse.h"  #include "parse.h"
   
 void Peval(), Psetprec(), Ptodouble();  void Peval(), Psetprec(), Psetbprec(), Ptodouble(), Psetround();
   void Pmpfr_ai();
   void Pmpfr_eint(), Pmpfr_erf(),Pmpfr_li2();
   void Pmpfr_zeta();
   void Pmpfr_j0(), Pmpfr_j1();
   void Pmpfr_y0(), Pmpfr_y1();
   void Pmpfr_gamma(), Pmpfr_lngamma(), Pmpfr_digamma();
   void Pmpfr_floor(), Pmpfr_round(), Pmpfr_ceil();
   
 struct ftab bf_tab[] = {  struct ftab bf_tab[] = {
         {"eval",Peval,-2},          {"eval",Peval,-2},
         {"setprec",Psetprec,-1},          {"setprec",Psetprec,-1},
           {"setbprec",Psetbprec,-1},
           {"setround",Psetround,-1},
         {"todouble",Ptodouble,1},          {"todouble",Ptodouble,1},
           {"mpfr_ai",Pmpfr_ai,-2},
           {"mpfr_zeta",Pmpfr_zeta,-2},
           {"mpfr_j0",Pmpfr_j0,-2},
           {"mpfr_j1",Pmpfr_j1,-2},
           {"mpfr_y0",Pmpfr_y0,-2},
           {"mpfr_y1",Pmpfr_y1,-2},
           {"mpfr_eint",Pmpfr_eint,-2},
           {"mpfr_erf",Pmpfr_erf,-2},
           {"mpfr_li2",Pmpfr_li2,-2},
           {"mpfr_gamma",Pmpfr_gamma,-2},
           {"mpfr_lngamma",Pmpfr_gamma,-2},
           {"mpfr_digamma",Pmpfr_gamma,-2},
           {"mpfr_floor",Pmpfr_floor,-2},
           {"mpfr_ceil",Pmpfr_ceil,-2},
           {"mpfr_round",Pmpfr_round,-2},
         {0,0,0},          {0,0,0},
 };  };
   
   int mpfr_roundmode = MPFR_RNDN;
   
 void Ptodouble(NODE arg,Num *rp)  void Ptodouble(NODE arg,Num *rp)
 {  {
         double r,i;          double r,i;
Line 42  void Ptodouble(NODE arg,Num *rp)
Line 68  void Ptodouble(NODE arg,Num *rp)
         }          }
 }  }
   
 void Peval(arg,rp)  void Peval(NODE arg,Obj *rp)
 NODE arg;  
 Obj *rp;  
 {  {
   int prec;    int prec;
   
         asir_assert(ARG0(arg),O_R,"eval");          asir_assert(ARG0(arg),O_R,"eval");
   if ( argc(arg) == 2 ) {    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;      if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN;
     else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX;      else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX;
   } else    } else
     prec = 0;      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)  void Psetprec(NODE arg,Obj *rp)
 {  {
         int p;          int p;
         Q q;          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 ) {          if ( arg ) {
                 asir_assert(ARG0(arg),O_N,"setprec");                  asir_assert(ARG0(arg),O_N,"setbprec");
                 p = QTOS((Q)ARG0(arg));                  prec = QTOS((Q)ARG0(arg));
                 if ( p > 0 )                  if ( p > 0 )
                         prec = p;                          prec = p;
         }          }
Line 79  void Psetprec(NODE arg,Obj *rp)
Line 126  void Psetprec(NODE arg,Obj *rp)
         mpfr_set_default_prec(prec);          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);  Num tobf(Num a,int prec);
   
 void mp_pi(NODE arg,BF *rp)  void mp_pi(NODE arg,BF *rp)
Line 89  void mp_pi(NODE arg,BF *rp)
Line 173  void mp_pi(NODE arg,BF *rp)
         prec = arg ? QTOS((Q)ARG0(arg)) : 0;          prec = arg ? QTOS((Q)ARG0(arg)) : 0;
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 103  void mp_e(NODE arg,BF *rp)
Line 187  void mp_e(NODE arg,BF *rp)
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
         mpfr_init(one);          mpfr_init(one);
         mpfr_set_ui(one,1,MPFR_RNDN);          mpfr_set_ui(one,1,mpfr_roundmode);
         mpfr_exp(r->body,one,MPFR_RNDN);          mpfr_exp(r->body,one,mpfr_roundmode);
     *rp = r;      *rp = r;
 }  }
   
Line 118  void mp_sin(NODE arg,BF *rp)
Line 202  void mp_sin(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 132  void mp_cos(NODE arg,BF *rp)
Line 216  void mp_cos(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 146  void mp_tan(NODE arg,BF *rp)
Line 230  void mp_tan(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 160  void mp_asin(NODE arg,BF *rp)
Line 244  void mp_asin(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
 void mp_acos(NODE arg,BF *rp)  void mp_acos(NODE arg,BF *rp)
Line 173  void mp_acos(NODE arg,BF *rp)
Line 257  void mp_acos(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
 void mp_atan(NODE arg,BF *rp)  void mp_atan(NODE arg,BF *rp)
Line 186  void mp_atan(NODE arg,BF *rp)
Line 270  void mp_atan(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 200  void mp_sinh(NODE arg,BF *rp)
Line 284  void mp_sinh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 214  void mp_cosh(NODE arg,BF *rp)
Line 298  void mp_cosh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 228  void mp_tanh(NODE arg,BF *rp)
Line 312  void mp_tanh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 242  void mp_asinh(NODE arg,BF *rp)
Line 326  void mp_asinh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
 void mp_acosh(NODE arg,BF *rp)  void mp_acosh(NODE arg,BF *rp)
Line 255  void mp_acosh(NODE arg,BF *rp)
Line 339  void mp_acosh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
 void mp_atanh(NODE arg,BF *rp)  void mp_atanh(NODE arg,BF *rp)
Line 268  void mp_atanh(NODE arg,BF *rp)
Line 352  void mp_atanh(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 282  void mp_exp(NODE arg,BF *rp)
Line 366  void mp_exp(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 296  void mp_log(NODE arg,BF *rp)
Line 380  void mp_log(NODE arg,BF *rp)
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
 }  }
   
Line 311  void mp_pow(NODE arg,BF *rp)
Line 395  void mp_pow(NODE arg,BF *rp)
         e = tobf(ARG1(arg),prec);          e = tobf(ARG1(arg),prec);
         NEWBF(r);          NEWBF(r);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);          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;      *rp = r;
   }
   
   #define SETPREC \
    (prec)=NEXT(arg)?QTOS((Q)ARG1(arg)):0;\
    (prec)*=3.32193;\
    (a)=tobf(ARG0(arg),prec);\
    NEWBF(r);\
    prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
   
   
   void Pmpfr_gamma(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_gamma(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_lngamma(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_lngamma(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_digamma(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_digamma(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_zeta(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_zeta(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_eint(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_eint(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_erf(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_erf(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_j0(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_j0(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_j1(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_j1(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_y0(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_y0(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_y1(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_y1(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_li2(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_li2(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_ai(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_ai(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_floor(NODE arg,Q *rp)
   {
           Num a;
     int prec;
           BF r;
           mpz_t t;
           GZ rz;
   
     SETPREC
           mpfr_floor(r->body,((BF)a)->body);
           mpz_init(t);
           mpfr_get_z(t,r->body,mpfr_roundmode);
           MPZTOGZ(t,rz);
           *rp = gztoz(rz);
   }
   
   void Pmpfr_ceil(NODE arg,Q *rp)
   {
           Num a;
     int prec;
           BF r;
           mpz_t t;
           GZ rz;
   
     SETPREC
           mpfr_ceil(r->body,((BF)a)->body);
           mpz_init(t);
           mpfr_get_z(t,r->body,mpfr_roundmode);
           MPZTOGZ(t,rz);
           *rp = gztoz(rz);
   }
   
   void Pmpfr_round(NODE arg,Q *rp)
   {
           Num a;
     int prec;
           BF r;
           mpz_t t;
           GZ rz;
   
     SETPREC
           mpfr_round(r->body,((BF)a)->body);
           mpz_init(t);
           mpfr_get_z(t,r->body,mpfr_roundmode);
           MPZTOGZ(t,rz);
           *rp = gztoz(rz);
 }  }
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.8
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>