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

Diff for /OpenXM_contrib2/asir2000/builtin/math.c between version 1.2 and 1.13

version 1.2, 2000/08/21 08:31:20 version 1.13, 2018/03/29 01:32:50
Line 23 
Line 23 
  * shall be made on your publication or presentation in any form of the   * shall be made on your publication or presentation in any form of the
  * results obtained by use of the SOFTWARE.   * results obtained by use of the SOFTWARE.
  * (4) In the event that you modify the SOFTWARE, you shall notify FLL by   * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
  * e-mail at risa-admin@flab.fujitsu.co.jp of the detailed specification   * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
  * for such modification or the source code of the modified part of the   * for such modification or the source code of the modified part of the
  * SOFTWARE.   * SOFTWARE.
  *   *
Line 45 
Line 45 
  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,   * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.   * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
  *   *
  * $OpenXM: OpenXM_contrib2/asir2000/builtin/math.c,v 1.1.1.1 1999/12/03 07:39:07 noro Exp $   * $OpenXM: OpenXM_contrib2/asir2000/builtin/math.c,v 1.12 2015/08/14 13:51:54 fujimoto Exp $
 */  */
 #include "ca.h"  #include "ca.h"
 #include <math.h>  #include <math.h>
 #include "parse.h"  #include "parse.h"
   #if defined(VISUAL) || defined(__MINGW32__)
   #include <float.h>
   #endif
   
 void Pdsqrt(),Pdsin(),Pdcos(),Pdtan(),Pdasin(),Pdacos(),Pdatan(),Pdlog(),Pdexp();  void Pdsqrt(),Pdsin(),Pdcos(),Pdtan(),Pdasin(),Pdacos(),Pdatan(),Pdlog(),Pdexp();
 void Pabs(),Pdfloor(),Pdceil(),Pdrint();  void Pabs(),Pdfloor(),Pdceil(),Pdrint(),Pdisnan();
   
 struct ftab math_tab[] = {  struct ftab math_tab[] = {
         {"dsqrt",Pdsqrt,1},    {"dsqrt",Pdsqrt,1},
         {"dabs",Pabs,1},    {"dabs",Pabs,1},
         {"dsin",Pdsin,1},    {"dsin",Pdsin,1},
         {"dcos",Pdcos,1},    {"dcos",Pdcos,1},
         {"dtan",Pdtan,1},    {"dtan",Pdtan,1},
         {"dlog",Pdlog,1},    {"dlog",Pdlog,1},
         {"dexp",Pdexp,1},    {"dexp",Pdexp,1},
         {"dasin",Pdasin,1},    {"dasin",Pdasin,1},
         {"dacos",Pdacos,1},    {"dacos",Pdacos,1},
         {"datan",Pdatan,1},    {"datan",Pdatan,1},
         {"dfloor",Pdfloor,1},    {"floor",Pdfloor,1},
         {"dceil",Pdceil,1},    {"dfloor",Pdfloor,1},
         {"drint",Pdrint,1},    {"ceil",Pdceil,1},
         {0,0,0},    {"dceil",Pdceil,1},
     {"rint",Pdrint,1},
     {"drint",Pdrint,1},
     {"disnan",Pdisnan,1},
     {0,0,0},
 };  };
   
   void get_ri(Num z,double *r,double *i)
   {
     if ( !z ) {
       *r = 0; *i = 0; return;
     }
     if ( OID(z) != O_N )
       error("get_ri : invalid argument");
     switch ( NID(z) ) {
       case N_Q: case N_R: case N_B:
         *r = ToReal(z); *i = 0;
         break;
       case N_C:
         *r = ToReal(((C)z)->r);
         *i = ToReal(((C)z)->i);
         break;
       default:
         error("get_ri : invalid argument");
         break;
     }
   }
   
 void Pabs(arg,rp)  void Pabs(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s,r,i;
   
         s = fabs(ToReal(ARG0(arg)));    if ( !ARG0(arg) ) {
         MKReal(s,*rp);      *rp = 0; return;
     }
     get_ri((Num)ARG0(arg),&r,&i);
     if ( i == 0 )
       s = fabs(r);
     else if ( r == 0 )
       s = fabs(i);
     else
       s = sqrt(r*r+i*i);
     MKReal(s,*rp);
 }  }
   
 void Pdsqrt(arg,rp)  void Pdsqrt(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Num *rp;
 {  {
         double s;    double s,r,i,a;
     C z;
     Real real;
   
         s = sqrt(ToReal(ARG0(arg)));    if ( !ARG0(arg) ) {
         MKReal(s,*rp);      *rp = 0; return;
     }
     get_ri((Num)ARG0(arg),&r,&i);
     if ( i == 0 )
       if ( r > 0 ) {
         s = sqrt(r);
         MKReal(s,real);
         *rp = (Num)real;
       } else {
         NEWC(z);
         z->r = 0;
         s = sqrt(-r); MKReal(s,real); z->i = (Num)real;
         *rp = (Num)z;
       }
     else {
       a = sqrt(r*r+i*i);
       NEWC(z);
       s = sqrt((r+a)/2); MKReal(s,real); z->r = (Num)real;
       s = i>0?sqrt((-r+a)/2):-sqrt((-r+a)/2);
       MKReal(s,real); z->i = (Num)real;
       *rp = (Num)z;
     }
 }  }
   
 void Pdsin(arg,rp)  void Pdsin(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = sin(ToReal(ARG0(arg)));    s = sin(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdcos(arg,rp)  void Pdcos(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = cos(ToReal(ARG0(arg)));    s = cos(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdtan(arg,rp)  void Pdtan(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = tan(ToReal(ARG0(arg)));    s = tan(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdasin(arg,rp)  void Pdasin(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = asin(ToReal(ARG0(arg)));    s = asin(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdacos(arg,rp)  void Pdacos(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = acos(ToReal(ARG0(arg)));    s = acos(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdatan(arg,rp)  void Pdatan(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = atan(ToReal(ARG0(arg)));    s = atan(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdlog(arg,rp)  void Pdlog(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = log(ToReal(ARG0(arg)));    s = log(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdexp(arg,rp)  void Pdexp(arg,rp)
 NODE arg;  NODE arg;
 Real *rp;  Real *rp;
 {  {
         double s;    double s;
   
         s = exp(ToReal(ARG0(arg)));    s = exp(ToReal(ARG0(arg)));
         MKReal(s,*rp);    MKReal(s,*rp);
 }  }
   
 void Pdfloor(arg,rp)  void Pdfloor(arg,rp)
 NODE arg;  NODE arg;
 Q *rp;  Q *rp;
 {  {
         L a;    L a;
         unsigned int au,al;    unsigned int au,al;
         int sgn;    int sgn;
         Q q;    Q q;
         double d;    double d;
   
         if ( !ARG0(arg) ) {    if ( !ARG0(arg) ) {
                 *rp = 0;      *rp = 0;
                 return;      return;
         }    }
         d = floor(ToReal(ARG0(arg)));    d = floor(ToReal(ARG0(arg)));
         if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )    if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                 error("dfloor : OverFlow");      error("dfloor : OverFlow");
         a = (L)d;    if ( !d ) {
         if ( a < 0 ) {      *rp = 0;
                 sgn = -1;      return;
                 a = -a;    }
         } else    a = (L)d;
                 sgn = 1;    if ( a < 0 ) {
 #if defined(i386) || defined(__alpha) || defined(VISUAL)      sgn = -1;
         au = ((unsigned int *)&a)[1];      a = -a;
         al = ((unsigned int *)&a)[0];    } else
       sgn = 1;
   #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__x86_64)
     au = ((unsigned int *)&a)[1];
     al = ((unsigned int *)&a)[0];
 #else  #else
         al = ((unsigned int *)&a)[1];    al = ((unsigned int *)&a)[1];
         au = ((unsigned int *)&a)[0];    au = ((unsigned int *)&a)[0];
 #endif  #endif
         if ( au ) {    if ( au ) {
                 NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;      NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                 PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;      PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
         } else {    } else {
                 UTOQ(al,q); SGN(q) = sgn;      UTOQ(al,q); SGN(q) = sgn;
         }    }
         *rp = q;    *rp = q;
 }  }
   
 void Pdceil(arg,rp)  void Pdceil(arg,rp)
 NODE arg;  NODE arg;
 Q *rp;  Q *rp;
 {  {
         L a;    L a;
         unsigned int au,al;    unsigned int au,al;
         int sgn;    int sgn;
         Q q;    Q q;
         double d;    double d;
   
         if ( !ARG0(arg) ) {    if ( !ARG0(arg) ) {
                 *rp = 0;      *rp = 0;
                 return;      return;
         }    }
         d = ceil(ToReal(ARG0(arg)));    d = ceil(ToReal(ARG0(arg)));
         if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )    if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                 error("dceil : OverFlow");      error("dceil : OverFlow");
         a = (L)d;    if ( !d ) {
         if ( a < 0 ) {      *rp = 0;
                 sgn = -1;      return;
                 a = -a;    }
         } else    a = (L)d;
                 sgn = 1;    if ( a < 0 ) {
 #if defined(i386) || defined(__alpha) || defined(VISUAL)      sgn = -1;
         au = ((unsigned int *)&a)[1];      a = -a;
         al = ((unsigned int *)&a)[0];    } else
       sgn = 1;
   #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__x86_64)
     au = ((unsigned int *)&a)[1];
     al = ((unsigned int *)&a)[0];
 #else  #else
         al = ((unsigned int *)&a)[1];    al = ((unsigned int *)&a)[1];
         au = ((unsigned int *)&a)[0];    au = ((unsigned int *)&a)[0];
 #endif  #endif
         if ( au ) {    if ( au ) {
                 NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;      NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                 PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;      PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
         } else {    } else {
                 UTOQ(al,q); SGN(q) = sgn;      UTOQ(al,q); SGN(q) = sgn;
         }    }
         *rp = q;    *rp = q;
 }  }
   
 void Pdrint(arg,rp)  void Pdrint(arg,rp)
 NODE arg;  NODE arg;
 Q *rp;  Q *rp;
 {  {
         L a;    L a;
         unsigned int au,al;    unsigned int au,al;
         int sgn;    int sgn;
         Q q;    Q q;
         double d;    double d;
   
         if ( !ARG0(arg) ) {    if ( !ARG0(arg) ) {
                 *rp = 0;      *rp = 0;
                 return;      return;
         }    }
 #if defined(VISUAL)  #if defined(VISUAL) || defined(__MINGW32__)
         d = ToReal(ARG0(arg));    d = ToReal(ARG0(arg));
         if ( d > 0 )    if ( d > 0 )
                 d = floor(d+0.5);      d = floor(d+0.5);
         else    else
                 d = ceil(d-0.5);      d = ceil(d-0.5);
 #else  #else
         d = rint(ToReal(ARG0(arg)));    d = rint(ToReal(ARG0(arg)));
 #endif  #endif
         if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )    if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                 error("drint : OverFlow");      error("drint : OverFlow");
         a = (L)d;    a = (L)d;
         if ( a < 0 ) {    if ( a < 0 ) {
                 sgn = -1;      sgn = -1;
                 a = -a;      a = -a;
         } else    } else
                 sgn = 1;      sgn = 1;
 #if defined(i386) || defined(__alpha) || defined(VISUAL)  #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__x86_64)
         au = ((unsigned int *)&a)[1];    au = ((unsigned int *)&a)[1];
         al = ((unsigned int *)&a)[0];    al = ((unsigned int *)&a)[0];
 #else  #else
         al = ((unsigned int *)&a)[1];    al = ((unsigned int *)&a)[1];
         au = ((unsigned int *)&a)[0];    au = ((unsigned int *)&a)[0];
 #endif  #endif
         if ( au ) {    if ( au ) {
                 NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;      NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                 PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;      PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
         } else {    } else if ( al ) {
                 UTOQ(al,q); SGN(q) = sgn;      UTOQ(al,q); SGN(q) = sgn;
         }    } else
         *rp = q;      q = 0;
     *rp = q;
   }
   
   void Pdisnan(NODE arg,Q *rp)
   {
     Real r;
     double d;
   #if defined(VISUAL) || defined(__MINGW32__)
       int c;
   #endif
   
     r = (Real)ARG0(arg);
     if ( !r || !NUM(r) || !REAL(r) ) {
       *rp = 0;
       return;
     }
     d = ToReal(r);
   #if defined(VISUAL) || defined(__MINGW32__)
     c = _fpclass(d);
     if ( c == _FPCLASS_SNAN || c == _FPCLASS_QNAN ) *rp = ONE;
     else if ( c == _FPCLASS_PINF || c == _FPCLASS_NINF ) STOQ(2,*rp);
   #else
     if ( isnan(d) ) *rp = ONE;
     else if ( isinf(d) ) STOQ(2,*rp);
   #endif
     else *rp = 0;
 }  }
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.13
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>