[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.8 and 1.11

version 1.8, 2015/08/16 03:12:09 version 1.11, 2015/08/25 18:41:03
Line 1 
Line 1 
 /* $OpenXM: OpenXM_contrib2/asir2000/builtin/bfaux.c,v 1.7 2015/08/07 08:00:30 takayama Exp $ */  /* $OpenXM: OpenXM_contrib2/asir2000/builtin/bfaux.c,v 1.10 2015/08/20 08:42:07 noro Exp $ */
 #include "ca.h"  #include "ca.h"
 #include "parse.h"  #include "parse.h"
   
 void Peval(), Psetprec(), Psetbprec(), Ptodouble(), Psetround();  void Peval(), Psetprec(), Psetbprec(), Ptodouble(), Psetround();
 void Pmpfr_ai();  void Pmpfr_ai();
 void Pmpfr_eint(), Pmpfr_erf(),Pmpfr_li2();  void Pmpfr_eint(), Pmpfr_erf(), Pmpfr_erfc(), Pmpfr_li2();
 void Pmpfr_zeta();  void Pmpfr_zeta();
 void Pmpfr_j0(), Pmpfr_j1();  void Pmpfr_j0(), Pmpfr_j1();
 void Pmpfr_y0(), Pmpfr_y1();  void Pmpfr_y0(), Pmpfr_y1();
Line 25  struct ftab bf_tab[] = {
Line 25  struct ftab bf_tab[] = {
         {"mpfr_y1",Pmpfr_y1,-2},          {"mpfr_y1",Pmpfr_y1,-2},
         {"mpfr_eint",Pmpfr_eint,-2},          {"mpfr_eint",Pmpfr_eint,-2},
         {"mpfr_erf",Pmpfr_erf,-2},          {"mpfr_erf",Pmpfr_erf,-2},
           {"mpfr_erfc",Pmpfr_erfc,-2},
         {"mpfr_li2",Pmpfr_li2,-2},          {"mpfr_li2",Pmpfr_li2,-2},
         {"mpfr_gamma",Pmpfr_gamma,-2},          {"mpfr_gamma",Pmpfr_gamma,-2},
         {"mpfr_lngamma",Pmpfr_gamma,-2},          {"mpfr_lngamma",Pmpfr_gamma,-2},
Line 96  void Psetprec(NODE arg,Obj *rp)
Line 97  void Psetprec(NODE arg,Obj *rp)
   STOQ(dprec,q); *rp = (Obj)q;    STOQ(dprec,q); *rp = (Obj)q;
         if ( arg ) {          if ( arg ) {
                 asir_assert(ARG0(arg),O_N,"setprec");                  asir_assert(ARG0(arg),O_N,"setprec");
                 prec = QTOS((Q)ARG0(arg))*3.32193;                  p = QTOS((Q)ARG0(arg))*3.32193;
                 if ( p > 0 )                  if ( p > 0 )
                         prec = p;                          prec = p;
         }          }
Line 117  void Psetbprec(NODE arg,Obj *rp)
Line 118  void Psetbprec(NODE arg,Obj *rp)
   STOQ(prec,q); *rp = (Obj)q;    STOQ(prec,q); *rp = (Obj)q;
         if ( arg ) {          if ( arg ) {
                 asir_assert(ARG0(arg),O_N,"setbprec");                  asir_assert(ARG0(arg),O_N,"setbprec");
                 prec = QTOS((Q)ARG0(arg));                  p = QTOS((Q)ARG0(arg));
                 if ( p > 0 )                  if ( p > 0 )
                         prec = p;                          prec = p;
         }          }
Line 167  Num tobf(Num a,int prec);
Line 168  Num tobf(Num a,int prec);
   
 void mp_pi(NODE arg,BF *rp)  void mp_pi(NODE arg,BF *rp)
 {  {
     int prec;    int prec;
         BF r;          BF r;
   
         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_roundmode);          mpfr_const_pi(r->body,mpfr_roundmode);
     *rp = r;    if ( !cmpbf((Num)r,0) ) r = 0;
     *rp = r;
 }  }
   
 void mp_e(NODE arg,BF *rp)  void mp_e(NODE arg,BF *rp)
 {  {
     int prec;    int prec;
         mpfr_t one;          mpfr_t one;
         BF r;          BF r;
   
Line 189  void mp_e(NODE arg,BF *rp)
Line 191  void mp_e(NODE arg,BF *rp)
         mpfr_init(one);          mpfr_init(one);
         mpfr_set_ui(one,1,mpfr_roundmode);          mpfr_set_ui(one,1,mpfr_roundmode);
         mpfr_exp(r->body,one,mpfr_roundmode);          mpfr_exp(r->body,one,mpfr_roundmode);
     *rp = r;    if ( !cmpbf((Num)r,0) ) r = 0;
     *rp = r;
 }  }
   
 void mp_sin(NODE arg,BF *rp)  void mpfr_or_mpc(NODE arg,int (*mpfr_f)(),int (*mpc_f)(),Num *rp)
 {  {
         Num a;          Num a;
     int prec;    int prec;
         BF r;          BF r,re,im;
     C c;
     mpc_t mpc,a1;
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;          prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : mpfr_get_default_prec();
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         NEWBF(r);    if ( a && NID(a)==N_C ) {
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);      mpc_init2(mpc,prec); mpc_init2(a1,prec);
         mpfr_sin(r->body,((BF)a)->body,mpfr_roundmode);      re = (BF)((C)a)->r; im = (BF)((C)a)->i;
     *rp = r;      mpc_set_fr_fr(a1,re->body,im->body,mpfr_roundmode);
       (*mpc_f)(mpc,a1,mpfr_roundmode);
       MPFRTOBF(mpc_realref(mpc),re);
       MPFRTOBF(mpc_imagref(mpc),im);
       if ( !cmpbf((Num)re,0) ) re = 0;
       if ( !cmpbf((Num)im,0) ) im = 0;
       if ( !im )
         *rp = (Num)re;
       else {
         NEWC(c); c->r = (Num)re; c->i = (Num)im;
         *rp = (Num)c;
       }
     } else {
             NEWBF(r);
             mpfr_init2(r->body,prec);
             (*mpfr_f)(r->body,((BF)a)->body,mpfr_roundmode);
       if ( !cmpbf((Num)r,0) ) r = 0;
       *rp = (Num)r;
     }
 }  }
   
 void mp_cos(NODE arg,BF *rp)  void mp_sin(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_sin,mpc_sin,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_tan(NODE arg,BF *rp)  void mp_cos(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_cos,mpc_cos,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_asin(NODE arg,BF *rp)  void mp_tan(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_tan,mpc_tan,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
 void mp_acos(NODE arg,BF *rp)  
 {  
         Num a;  
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  void mp_asin(NODE arg,Num *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_roundmode);  
     *rp = r;  
 }  
 void mp_atan(NODE arg,BF *rp)  
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_asin,mpc_asin,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_sinh(NODE arg,BF *rp)  void mp_acos(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_acos,mpc_acos,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_cosh(NODE arg,BF *rp)  void mp_atan(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_atan,mpc_atan,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_tanh(NODE arg,BF *rp)  void mp_sinh(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_sinh,mpc_sinh,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
   
 void mp_asinh(NODE arg,BF *rp)  void mp_cosh(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_cosh,mpc_cosh,rp);
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  
         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_roundmode);  
     *rp = r;  
 }  }
 void mp_acosh(NODE arg,BF *rp)  
 {  
         Num a;  
     int prec;  
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  void mp_tanh(NODE arg,Num *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_roundmode);  
     *rp = r;  
 }  
 void mp_atanh(NODE arg,BF *rp)  
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_tanh,mpc_tanh,rp);
     int prec;  }
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  void mp_asinh(NODE arg,Num *rp)
         a = tobf(ARG0(arg),prec);  {
         NEWBF(r);    mpfr_or_mpc(arg,mpfr_asinh,mpc_asinh,rp);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);  
         mpfr_atanh(r->body,((BF)a)->body,mpfr_roundmode);  
     *rp = r;  
 }  }
   
 void mp_exp(NODE arg,BF *rp)  void mp_acosh(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_acosh,mpc_acosh,rp);
     int prec;  }
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  void mp_atanh(NODE arg,Num *rp)
         a = tobf(ARG0(arg),prec);  {
         NEWBF(r);    mpfr_or_mpc(arg,mpfr_atanh,mpc_atanh,rp);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);  
         mpfr_exp(r->body,((BF)a)->body,mpfr_roundmode);  
     *rp = r;  
 }  }
   
 void mp_log(NODE arg,BF *rp)  void mp_exp(NODE arg,Num *rp)
 {  {
         Num a;    mpfr_or_mpc(arg,mpfr_exp,mpc_exp,rp);
     int prec;  }
         BF r;  
   
         prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;  void mp_log(NODE arg,Num *rp)
         a = tobf(ARG0(arg),prec);  {
         NEWBF(r);    mpfr_or_mpc(arg,mpfr_log,mpc_log,rp);
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);  
         mpfr_log(r->body,((BF)a)->body,mpfr_roundmode);  
     *rp = r;  
 }  }
   
 void mp_pow(NODE arg,BF *rp)  void mp_pow(NODE arg,Num *rp)
 {  {
         Num a,e;          Num a,e;
     int prec;    int prec;
         BF r;          BF r,re,im;
     C c;
     mpc_t mpc,a1,e1;
   
         prec = NEXT(NEXT(arg)) ? QTOS((Q)ARG2(arg)) : 0;          prec = NEXT(NEXT(arg)) ? QTOS((Q)ARG2(arg)) : mpfr_get_default_prec();
         a = tobf(ARG0(arg),prec);          a = tobf(ARG0(arg),prec);
         e = tobf(ARG1(arg),prec);          e = tobf(ARG1(arg),prec);
         NEWBF(r);    if ( NID(a) == N_C || NID(e) == N_C || MPFR_SIGN(((BF)a)->body) < 0 ) {
         prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);      mpc_init2(mpc,prec); mpc_init2(a1,prec); mpc_init2(e1,prec);
         mpfr_pow(r->body,((BF)a)->body,((BF)e)->body,mpfr_roundmode);      if ( NID(a) == N_C ) {
     *rp = r;        re = (BF)((C)a)->r; im = (BF)((C)a)->i;
         mpc_set_fr_fr(a1,re->body,im->body,mpfr_roundmode);
       } else {
         re = (BF)a;
         mpc_set_fr(a1,re->body,mpfr_roundmode);
       }
       if ( NID(e) == N_C ) {
         re = (BF)((C)e)->r; im = (BF)((C)e)->i;
         mpc_set_fr_fr(e1,re->body,im->body,mpfr_roundmode);
       } else {
         re = (BF)e;
         mpc_set_fr(e1,re->body,mpfr_roundmode);
       }
       mpc_pow(mpc,a1,e1,mpfr_roundmode);
       MPFRTOBF(mpc_realref(mpc),re);
       MPFRTOBF(mpc_imagref(mpc),im);
       if ( !cmpbf((Num)re,0) ) re = 0;
       if ( !cmpbf((Num)im,0) ) im = 0;
       if ( !im )
         *rp = (Num)re;
       else {
         NEWC(c); c->r = (Num)re; c->i = (Num)im;
         *rp = (Num)c;
       }
     } else {
             NEWBF(r);
             mpfr_init2(r->body,prec);
             mpfr_pow(r->body,((BF)a)->body,((BF)e)->body,mpfr_roundmode);
       *rp = (Num)r;
     }
 }  }
   
 #define SETPREC \  #define SETPREC \
Line 470  void Pmpfr_erf(NODE arg,BF *rp)
Line 416  void Pmpfr_erf(NODE arg,BF *rp)
   
   SETPREC    SETPREC
         mpfr_erf(r->body,((BF)a)->body,mpfr_roundmode);          mpfr_erf(r->body,((BF)a)->body,mpfr_roundmode);
     *rp = r;
   }
   
   void Pmpfr_erfc(NODE arg,BF *rp)
   {
           Num a;
     int prec;
           BF r;
   
     SETPREC
           mpfr_erfc(r->body,((BF)a)->body,mpfr_roundmode);
   *rp = r;    *rp = r;
 }  }
   
Legend:
Removed from v.1.8  
changed lines
  Added in v.1.11
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>