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

Diff for /OpenXM_contrib2/asir2000/engine/f-itv.c between version 1.1 and 1.7

version 1.1, 2000/12/22 10:03:28 version 1.7, 2009/03/27 14:42:29
Line 1 
Line 1 
 /*  /*
  * $OpenXM: $   * $OpenXM: OpenXM_contrib2/asir2000/engine/f-itv.c,v 1.6 2005/01/11 07:12:51 saito Exp $
 */  */
 #if defined(INTERVAL)  #if defined(INTERVAL)
 #include "ca.h"  #include "ca.h"
 #include "base.h"  #include "base.h"
 #if PARI  #if defined(PARI)
 #include "genpari.h"  #include "genpari.h"
 #include "itv-pari.h"  #include "itv-pari.h"
 extern long prec;  long get_pariprec();
 #endif  #endif
   
 void ToBf(Num a, BF *rp)  void ToBf(Num a, BF *rp)
Line 24  void ToBf(Num a, BF *rp)
Line 24  void ToBf(Num a, BF *rp)
                 case N_Q:                  case N_Q:
                         ltop = avma;                          ltop = avma;
                         ritopa_i(NM((Q)a), SGN((Q)a), &pa);                          ritopa_i(NM((Q)a), SGN((Q)a), &pa);
                         pb = cgetr(prec);                          pb = cgetr(get_pariprec());
                         mpaff(pa, pb);                          mpaff(pa, pb);
                         if ( INT((Q)a) ) {                          if ( INT((Q)a) ) {
                                 lbot = avma;                                  lbot = avma;
Line 35  void ToBf(Num a, BF *rp)
Line 35  void ToBf(Num a, BF *rp)
                         } else {                          } else {
                                 patori(pb, &bn);                                  patori(pb, &bn);
                                 ritopa_i(DN((Q)a), 1, &pa);                                  ritopa_i(DN((Q)a), 1, &pa);
                                 pb = cgetr(prec);                                  pb = cgetr(get_pariprec());
                                 mpaff(pa, pb);                                  mpaff(pa, pb);
                                 lbot = avma;                                  lbot = avma;
                                 pb = gerepile(ltop, lbot, pb);                                  pb = gerepile(ltop, lbot, pb);
Line 94  void getulp(BF a, BF *au)
Line 94  void getulp(BF a, BF *au)
         cgiv(b);          cgiv(b);
 }  }
   
 void addulp(ItvF a, ItvF *c)  void addulp(IntervalBigFloat a, IntervalBigFloat *c)
 {  {
         Num     ai, as, aiu, asu, inf, sup;          Num     ai, as, aiu, asu, inf, sup;
   
Line 110  void addulp(ItvF a, ItvF *c)
Line 110  void addulp(ItvF a, ItvF *c)
         istoitv(inf,sup, (Itv *)c);          istoitv(inf,sup, (Itv *)c);
 }  }
   
 void additvf(ItvF a, ItvF b, ItvF *c)  void additvf(IntervalBigFloat a, IntervalBigFloat b, IntervalBigFloat *c)
 {  {
         Num     ai,as,bi,bs,mas,mbs,tmp;          Num     ai,as,bi,bs,mas,mbs,tmp;
         Num     inf,sup;          Num     inf,sup;
Line 135  printexpr(CO, bi);
Line 135  printexpr(CO, bi);
 printexpr(CO, bs);  printexpr(CO, bs);
 #endif  #endif
   
 #if 1  
                 addnum(0,ai,bi,&inf);                  addnum(0,ai,bi,&inf);
                 addnum(0,as,bs,&sup);                  addnum(0,as,bs,&sup);
                 istoitv(inf,sup,(Itv *)&tmp);                  istoitv(inf,sup,(Itv *)&tmp);
                 addulp((ItvF)tmp, c);                  addulp((IntervalBigFloat)tmp, c);
                 return;                  return;
 #else  
                 ltop = avma;  
                 ritopa(ai,&pa);  
                 ritopa(bi,&pb);  
                 lbot = avma;  
                 z = gerepile(ltop,lbot,PariAddDown(pa,pb));  
                 patori(z,&inf); cgiv(z);  
   
         /* MUST check if ai, as, bi, bs are bigfloat. */  
   
                 /*  as + bs = ( - ( (-as) + (-bs) ) ) */  
                 chsgnbf(as,&mas);  
                 chsgnbf(bs,&mbs);  
                 ltop = avma;  
                 ritopa(mas,&pa);  
                 ritopa(mbs,&pb);  
                 lbot = avma;  
                 z = gerepile(ltop,lbot,PariAddDown(pa,pb));  
                 patori(z,&tmp); cgiv(z);  
   
                 chsgnbf(tmp,&sup);  
                 istoitv(inf,sup,c);  
 #endif  
         }          }
 }  }
   
 void subitvf(ItvF a, ItvF b, ItvF *c)  void subitvf(IntervalBigFloat a, IntervalBigFloat b, IntervalBigFloat *c)
 {  {
         Num     ai,as,bi,bs,mas, mbs;          Num     ai,as,bi,bs,mas, mbs;
         Num     inf,sup,tmp;          Num     inf,sup,tmp;
Line 185  void subitvf(ItvF a, ItvF b, ItvF *c)
Line 161  void subitvf(ItvF a, ItvF b, ItvF *c)
                 ToBf(inf, (BF *)&ai); ToBf(sup, (BF *)&as);                  ToBf(inf, (BF *)&ai); ToBf(sup, (BF *)&as);
                 itvtois((Itv)b,&inf,&sup);                  itvtois((Itv)b,&inf,&sup);
                 ToBf(inf, (BF *)&bi); ToBf(sup, (BF *)&bs);                  ToBf(inf, (BF *)&bi); ToBf(sup, (BF *)&bs);
 #if 1  
                 subnum(0,ai,bs,&inf);                  subnum(0,ai,bs,&inf);
                 subnum(0,as,bi,&sup);                  subnum(0,as,bi,&sup);
                 istoitv(inf,sup,(Itv *)&tmp);                  istoitv(inf,sup,(Itv *)&tmp);
                 addulp((ItvF)tmp, c);                  addulp((IntervalBigFloat)tmp, c);
 #else  
   
 /* MUST check if ai, as, bi, bs are bigfloat. */  
                 /* ai - bs = ai + (-bs)  */  
                 chsgnbf(bs,&mbs);  
                 ltop = avma;  
                 ritopa(ai,&pa);  
                 ritopa(mbs,&pb);  
                 lbot = avma;  
                 z = gerepile(ltop,lbot,PariAddDown(pa,pb));  
                 patori(z,&inf); cgiv(z);  
   
                 /* as - bi = ( - ( bi + (-as) ) ) */  
                 chsgnbf(as,&mas);  
                 ltop = avma;  
                 ritopa(mas,&pa);  
                 ritopa(bi,&pb);  
                 lbot = avma;  
                 z = gerepile(ltop,lbot,PariAddDown(pa,pb));  
                 patori(z,&tmp); cgiv(z);  
   
                 chsgnbf(tmp,&sup);  
                 istoitv(inf,sup,c);  
 #endif  
         }          }
 }  }
   
 void mulitvf(ItvF a, ItvF b, ItvF *c)  void mulitvf(IntervalBigFloat a, IntervalBigFloat b, IntervalBigFloat *c)
 {  {
         Num     ai,as,bi,bs,a1,a2,b1,b2,c1,c2,p,t,tmp;          Num     ai,as,bi,bs,a1,a2,b1,b2,c1,c2,p,t,tmp;
         Num     inf, sup;          Num     inf, sup;
Line 278  void mulitvf(ItvF a, ItvF b, ItvF *c)
Line 229  void mulitvf(ItvF a, ItvF b, ItvF *c)
                         subnum(0,0,c1,&t);                          subnum(0,0,c1,&t);
                         istoitv(c2,t,(Itv *)&tmp);                          istoitv(c2,t,(Itv *)&tmp);
                 }                  }
                 addulp((ItvF)tmp, c);                  addulp((IntervalBigFloat)tmp, c);
         }          }
 }  }
   
Line 302  int     itvinitvf(Itv a, Itv b)
Line 253  int     itvinitvf(Itv a, Itv b)
         else return 0;          else return 0;
 }  }
   
 void divitvf(ItvF a, ItvF b, ItvF *c)  void divitvf(IntervalBigFloat a, IntervalBigFloat b, IntervalBigFloat *c)
 {  {
         Num     ai,as,bi,bs,a1,a2,b1,b2,c1,c2,t,tmp;          Num     ai,as,bi,bs,a1,a2,b1,b2,c1,c2,t,tmp;
         Num     inf, sup;          Num     inf, sup;
Line 352  void divitvf(ItvF a, ItvF b, ItvF *c)
Line 303  void divitvf(ItvF a, ItvF b, ItvF *c)
                         subnum(0,0,c1,&t);                          subnum(0,0,c1,&t);
                         istoitv(c2,t,(Itv *)&tmp);                          istoitv(c2,t,(Itv *)&tmp);
                 }                  }
                 addulp((ItvF)tmp, c);                  addulp((IntervalBigFloat)tmp, c);
         }          }
 }  }
   
Line 368  void pwritvf(Itv a, Num e, Itv *c)
Line 319  void pwritvf(Itv a, Num e, Itv *c)
         else if ( NID(a) <= N_B )          else if ( NID(a) <= N_B )
                 pwrnum(0,(Num)a,e,(Num *)c);                  pwrnum(0,(Num)a,e,(Num *)c);
         else if ( !INT(e) ) {          else if ( !INT(e) ) {
 #if PARI && 0  #if defined(PARI) && 0
                 GEN pa,pe,z;                  gpui_ri((Obj)a,(Obj)c,(Obj *)c);
                 int ltop,lbot;  
   
                 ltop = avma; ritopa(a,&pa); ritopa(e,&pe); lbot = avma;  
                 z = gerepile(ltop,lbot,gpui(pa,pe,prec));  
                 patori(z,c); cgiv(z);  
 #else  #else
                 error("pwritv : can't calculate a fractional power");                  error("pwritv : can't calculate a fractional power");
 #endif  #endif
Line 392  void pwritv0f(Itv a, int e, Itv *c)
Line 338  void pwritv0f(Itv a, int e, Itv *c)
 {  {
         Num inf, sup;          Num inf, sup;
         Num ai,Xmin,Xmax;          Num ai,Xmin,Xmax;
         ItvF tmp;          IntervalBigFloat tmp;
         Q       ne;          Q       ne;
   
         if ( e == 1 )          if ( e == 1 )
Line 426  void pwritv0f(Itv a, int e, Itv *c)
Line 372  void pwritv0f(Itv a, int e, Itv *c)
                 if ( ! Xmax )   sup = 0;                  if ( ! Xmax )   sup = 0;
                 else            pwrbf(Xmax,(Num)ne,&sup);                  else            pwrbf(Xmax,(Num)ne,&sup);
                 istoitv(inf,sup,(Itv *)&tmp);                  istoitv(inf,sup,(Itv *)&tmp);
                 addulp(tmp, (ItvF *)c);                  addulp(tmp, (IntervalBigFloat *)c);
         }          }
 }  }
   
Line 451  int cmpitvf(Itv a, Itv b)
Line 397  int cmpitvf(Itv a, Itv b)
         int     s,t;          int     s,t;
   
         if ( !a ) {          if ( !a ) {
                 if ( !b || (NID(b)<=N_B) )                  if ( !b || (NID(b)<=N_B) ) {
                         return compnum(0,(Num)a,(Num)b);                          return compnum(0,(Num)a,(Num)b);
                 else                  } else {
                         return -1;                          itvtois(b,&bi,&bs);
                           if ( compnum(0,(Num)a,bs) > 0 ) return 1;
                           else if ( compnum(0,bi,(Num)a) > 0 ) return -1;
                           else  return 0;
                   }
         } else if ( !b ) {          } else if ( !b ) {
                 if ( !a || (NID(a)<=N_B) )                  if ( !a || (NID(a)<=N_B) ) {
                         return compnum(0,(Num)a,(Num)b);                          return compnum(0,(Num)a,(Num)b);
                 else                  } else {
                         return initvp((Num)b,a);                          itvtois(a,&ai,&as);
                           if ( compnum(0,ai,(Num)b) > 0 ) return 1;
                           else if ( compnum(0,(Num)b,as) > 0 ) return -1;
                           else  return 0;
                   }
         } else {          } else {
                 itvtois(a,&ai,&as);                  itvtois(a,&ai,&as);
                 itvtois(b,&bi,&bs);                  itvtois(b,&bi,&bs);
                 s = compnum(0,ai,bi);                  s = compnum(0,ai,bs) ;
                 t = compnum(0,as,bs);                  t = compnum(0,bi,as) ;
                 if ( !s && !t ) return 0;                  if ( s > 0 ) return 1;
                 else  return -1;                  else if ( t > 0 ) return -1;
                   else  return 0;
         }          }
 }  }
   
Line 474  void miditvf(Itv a, Num *b)
Line 429  void miditvf(Itv a, Num *b)
 {  {
         Num     ai,as;          Num     ai,as;
         Num     t;          Num     t;
         Q       TWO;  
   
         if ( ! a ) *b = 0;          if ( ! a ) *b = 0;
         else if ( (NID(a) <= N_B) )          else if ( (NID(a) <= N_B) )
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.7
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>