[BACK]Return to d-itv.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2018 / engine

Diff for /OpenXM_contrib2/asir2018/engine/d-itv.c between version 1.3 and 1.4

version 1.3, 2019/06/04 07:11:23 version 1.4, 2019/10/17 03:03:12
Line 1 
Line 1 
 /*  /*
  * $OpenXM: OpenXM_contrib2/asir2018/engine/d-itv.c,v 1.2 2018/09/28 08:20:28 noro Exp $   * $OpenXM: OpenXM_contrib2/asir2018/engine/d-itv.c,v 1.3 2019/06/04 07:11:23 kondoh Exp $
 */  */
 #if defined(INTERVAL)  #if defined(INTERVAL)
 #include <float.h>  #include <float.h>
Line 11 
Line 11 
 #endif  #endif
 #endif  #endif
   
   double  addulpd(double);
   double  subulpd(double);
   extern int mpfr_roundmode;
   Num tobf(Num,int);
   
   
 #if defined(ITVDEBUG)  #if defined(ITVDEBUG)
 void printbinint(int d)  void printbinint(int d)
 {  {
Line 33  void printbinint(int d)
Line 39  void printbinint(int d)
 }  }
 #endif  #endif
   
   #if 0
 double NatToRealUp(N a, int *expo)  double NatToRealUp(N a, int *expo)
 {  {
   int *p;    int *p;
Line 181  static double  Q2doubleUp(Q a)
Line 188  static double  Q2doubleUp(Q a)
   }    }
 }  }
   
 #if 0  
 static double  PARI2doubleDown(BF a)  static double  PARI2doubleDown(BF a)
 {  {
         //GEN     p;          //GEN     p;
Line 200  static double  PARI2doubleUp(BF a)
Line 206  static double  PARI2doubleUp(BF a)
 }  }
 #endif  #endif
   
   static double  Q2doubleDown(Q a)
   {
           Num b;
           double rd;
           int current_roundmode=0;
   
           current_roundmode = mpfr_roundmode;
           mpfr_roundmode = MPFR_RNDD;
           b = tobf((Num)a,60);
           mpfr_roundmode = current_roundmode;
           rd = mpfr_get_d(BDY((BF)b),MPFR_RNDD);
           return rd;
   }
   
   static double  Q2doubleUp(Q a)
   {
           Num b;
           double rd;
           int current_roundmode=0;
   
           current_roundmode = mpfr_roundmode;
           mpfr_roundmode = MPFR_RNDU;
           b = tobf((Num)a,60);
           mpfr_roundmode = current_roundmode;
           rd = mpfr_get_d(BDY((BF)b),MPFR_RNDU);
           return rd;
   }
   
 double  subulpd(double d)  double  subulpd(double d)
 {  {
   double inf;    double inf;
Line 232  double  ToRealDown(Num a)
Line 266  double  ToRealDown(Num a)
     case N_R:      case N_R:
       inf = subulpd(BDY((Real)a)); break;        inf = subulpd(BDY((Real)a)); break;
     case N_B:      case N_B:
       //inf = PARI2doubleDown((BF)a); break;                  inf = mpfr_get_d(BDY((BF)a),MPFR_RNDD);
         inf = 0;                  break;
       error("ToRealDown: not supported operands.");  
         break;  
     case N_IP:      case N_IP:
       inf = ToRealDown(INF((Itv)a));        inf = ToRealDown(INF((Itv)a));
       break;        break;
Line 261  double  ToRealUp(Num a)
Line 293  double  ToRealUp(Num a)
     case N_R:      case N_R:
       sup = addulpd(BDY((Real)a)); break;        sup = addulpd(BDY((Real)a)); break;
     case N_B:      case N_B:
       //sup = PARI2doubleUp((BF)a); break;                  sup = mpfr_get_d(BDY((BF)a),MPFR_RNDU);
       sup = 0;                  break;
       error("ToRealUp: not supported operands.");  
       break;  
     case N_IP:      case N_IP:
       sup = ToRealUp(SUP((Itv)a)); break;        sup = ToRealUp(SUP((Itv)a)); break;
     case N_IntervalDouble:      case N_IntervalDouble:
Line 291  void  Num2double(Num a, double *inf, double *sup)
Line 321  void  Num2double(Num a, double *inf, double *sup)
       *sup = BDY((Real)a);        *sup = BDY((Real)a);
       break;        break;
     case N_B:      case N_B:
       //*inf = PARI2doubleDown((BF)a);  
       //*sup = PARI2doubleUp((BF)a);  
       *inf = mpfr_get_d(BDY((BF)a), MPFR_RNDD);        *inf = mpfr_get_d(BDY((BF)a), MPFR_RNDD);
       *sup = mpfr_get_d(BDY((BF)a), MPFR_RNDU);        *sup = mpfr_get_d(BDY((BF)a), MPFR_RNDU);
       break;        break;
Line 493  void  divitvd(Num a, Num b, IntervalDouble *c)
Line 521  void  divitvd(Num a, Num b, IntervalDouble *c)
   
 void  pwritvd(Num a, Num e, IntervalDouble *c)  void  pwritvd(Num a, Num e, IntervalDouble *c)
 {  {
   int  ei;    long  ei;
   IntervalDouble  t;    IntervalDouble  t;
   
   if ( !e ) {    if ( !e ) {
Line 511  void  pwritvd(Num a, Num e, IntervalDouble *c)
Line 539  void  pwritvd(Num a, Num e, IntervalDouble *c)
     error("pwritvd : can't calculate a fractional power");      error("pwritvd : can't calculate a fractional power");
 #endif  #endif
   } else {    } else {
     ei = QTOS((Q)e);      //ei = QTOS((Q)e);
           ei = mpz_get_si(BDY((Z)e));
     if (ei<0) ei = -ei;      if (ei<0) ei = -ei;
     pwritv0d((IntervalDouble)a,ei,&t);      pwritv0d((IntervalDouble)a,ei,&t);
     if ( SGN((Q)e) < 0 )  //    if ( SGN((Q)e) < 0 )
       if ( sgnq((Q)e) < 0 )
       divnum(0,(Num)ONE,(Num)t,(Num *)c);        divnum(0,(Num)ONE,(Num)t,(Num *)c);
     else      else
       *c = t;        *c = t;
   }    }
 }  }
   
 void  pwritv0d(IntervalDouble a, int e, IntervalDouble *c)  void  pwritv0d(IntervalDouble a, long e, IntervalDouble *c)
 {  {
   double inf, sup;    double inf, sup;
   double t, Xmin, Xmax;    double t, Xmin, Xmax;
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.4
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>