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

Annotation of OpenXM_contrib2/asir2018/builtin/math.c, Revision 1.2

1.1       noro        1: /*
                      2:  * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
                      3:  * All rights reserved.
                      4:  *
                      5:  * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
                      6:  * non-exclusive and royalty-free license to use, copy, modify and
                      7:  * redistribute, solely for non-commercial and non-profit purposes, the
                      8:  * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
                      9:  * conditions of this Agreement. For the avoidance of doubt, you acquire
                     10:  * only a limited right to use the SOFTWARE hereunder, and FLL or any
                     11:  * third party developer retains all rights, including but not limited to
                     12:  * copyrights, in and to the SOFTWARE.
                     13:  *
                     14:  * (1) FLL does not grant you a license in any way for commercial
                     15:  * purposes. You may use the SOFTWARE only for non-commercial and
                     16:  * non-profit purposes only, such as academic, research and internal
                     17:  * business use.
                     18:  * (2) The SOFTWARE is protected by the Copyright Law of Japan and
                     19:  * international copyright treaties. If you make copies of the SOFTWARE,
                     20:  * with or without modification, as permitted hereunder, you shall affix
                     21:  * to all such copies of the SOFTWARE the above copyright notice.
                     22:  * (3) An explicit reference to this SOFTWARE and its copyright owner
                     23:  * shall be made on your publication or presentation in any form of the
                     24:  * results obtained by use of the SOFTWARE.
                     25:  * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
                     26:  * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
                     27:  * for such modification or the source code of the modified part of the
                     28:  * SOFTWARE.
                     29:  *
                     30:  * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
                     31:  * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
                     32:  * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
                     33:  * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
                     34:  * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
                     35:  * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
                     36:  * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
                     37:  * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
                     38:  * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
                     39:  * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
                     40:  * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
                     41:  * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
                     42:  * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
                     43:  * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
                     44:  * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
                     45:  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
                     46:  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
                     47:  *
1.2     ! noro       48:  * $OpenXM: OpenXM_contrib2/asir2018/builtin/math.c,v 1.1 2018/09/19 05:45:06 noro Exp $
1.1       noro       49: */
                     50: #include "ca.h"
                     51: #include <math.h>
                     52: #include "parse.h"
                     53: #if defined(VISUAL) || defined(__MINGW32__)
                     54: #include <float.h>
                     55: #endif
                     56:
                     57: void get_ri(Num z,double *r,double *i);
                     58: void Pabs(NODE arg,Real *rp);
                     59: void Pdsqrt(NODE arg,Num *rp);
                     60: void Pdsin(NODE arg,Real *rp);
                     61: void Pdcos(NODE arg,Real *rp);
                     62: void Pdtan(NODE arg,Real *rp);
                     63: void Pdasin(NODE arg,Real *rp);
                     64: void Pdacos(NODE arg,Real *rp);
                     65: void Pdatan(NODE arg,Real *rp);
                     66: void Pdlog(NODE arg,Real *rp);
                     67: void Pdexp(NODE arg,Real *rp);
                     68: void Pdfloor(NODE arg,Z *rp);
                     69: void Pdceil(NODE arg,Z *rp);
                     70: void Pdrint(NODE arg,Z *rp);
                     71: void Pdisnan(NODE arg,Z *rp);
                     72: struct ftab math_tab[] = {
                     73:   {"dsqrt",Pdsqrt,1},
                     74:   {"dabs",Pabs,1},
                     75:   {"dsin",Pdsin,1},
                     76:   {"dcos",Pdcos,1},
                     77:   {"dtan",Pdtan,1},
                     78:   {"dlog",Pdlog,1},
                     79:   {"dexp",Pdexp,1},
                     80:   {"dasin",Pdasin,1},
                     81:   {"dacos",Pdacos,1},
                     82:   {"datan",Pdatan,1},
                     83:   {"floor",Pdfloor,1},
                     84:   {"dfloor",Pdfloor,1},
                     85:   {"ceil",Pdceil,1},
                     86:   {"dceil",Pdceil,1},
                     87:   {"rint",Pdrint,1},
                     88:   {"drint",Pdrint,1},
                     89:   {"disnan",Pdisnan,1},
                     90:   {0,0,0},
                     91: };
                     92:
                     93: void get_ri(Num z,double *r,double *i)
                     94: {
                     95:   if ( !z ) {
                     96:     *r = 0; *i = 0; return;
                     97:   }
                     98:   if ( OID(z) != O_N )
                     99:     error("get_ri : invalid argument");
                    100:   switch ( NID(z) ) {
                    101:     case N_Q: case N_R: case N_B:
                    102:       *r = ToReal(z); *i = 0;
                    103:       break;
                    104:     case N_C:
                    105:       *r = ToReal(((C)z)->r);
                    106:       *i = ToReal(((C)z)->i);
                    107:       break;
                    108:     default:
                    109:       error("get_ri : invalid argument");
                    110:       break;
                    111:   }
                    112: }
                    113:
                    114: void Pabs(NODE arg,Real *rp)
                    115: {
                    116:   double s,r,i;
                    117:
                    118:   if ( !ARG0(arg) ) {
                    119:     *rp = 0; return;
                    120:   }
                    121:   get_ri((Num)ARG0(arg),&r,&i);
                    122:   if ( i == 0 )
                    123:     s = fabs(r);
                    124:   else if ( r == 0 )
                    125:     s = fabs(i);
                    126:   else
                    127:     s = sqrt(r*r+i*i);
                    128:   MKReal(s,*rp);
                    129: }
                    130:
                    131: void Pdsqrt(NODE arg,Num *rp)
                    132: {
                    133:   double s,r,i,a;
                    134:   C z;
                    135:   Real real;
                    136:
                    137:   if ( !ARG0(arg) ) {
                    138:     *rp = 0; return;
                    139:   }
                    140:   get_ri((Num)ARG0(arg),&r,&i);
                    141:   if ( i == 0 )
                    142:     if ( r > 0 ) {
                    143:       s = sqrt(r);
                    144:       MKReal(s,real);
                    145:       *rp = (Num)real;
                    146:     } else {
                    147:       NEWC(z);
                    148:       z->r = 0;
                    149:       s = sqrt(-r); MKReal(s,real); z->i = (Num)real;
                    150:       *rp = (Num)z;
                    151:     }
                    152:   else {
                    153:     a = sqrt(r*r+i*i);
                    154:     NEWC(z);
                    155:     s = sqrt((r+a)/2); MKReal(s,real); z->r = (Num)real;
                    156:     s = i>0?sqrt((-r+a)/2):-sqrt((-r+a)/2);
                    157:     MKReal(s,real); z->i = (Num)real;
                    158:     *rp = (Num)z;
                    159:   }
                    160: }
                    161:
                    162: void Pdsin(NODE arg,Real *rp)
                    163: {
                    164:   double s;
                    165:
                    166:   s = sin(ToReal(ARG0(arg)));
                    167:   MKReal(s,*rp);
                    168: }
                    169:
                    170: void Pdcos(NODE arg,Real *rp)
                    171: {
                    172:   double s;
                    173:
                    174:   s = cos(ToReal(ARG0(arg)));
                    175:   MKReal(s,*rp);
                    176: }
                    177:
                    178: void Pdtan(NODE arg,Real *rp)
                    179: {
                    180:   double s;
                    181:
                    182:   s = tan(ToReal(ARG0(arg)));
                    183:   MKReal(s,*rp);
                    184: }
                    185:
                    186: void Pdasin(NODE arg,Real *rp)
                    187: {
                    188:   double s;
                    189:
                    190:   s = asin(ToReal(ARG0(arg)));
                    191:   MKReal(s,*rp);
                    192: }
                    193:
                    194: void Pdacos(NODE arg,Real *rp)
                    195: {
                    196:   double s;
                    197:
                    198:   s = acos(ToReal(ARG0(arg)));
                    199:   MKReal(s,*rp);
                    200: }
                    201:
                    202: void Pdatan(NODE arg,Real *rp)
                    203: {
                    204:   double s;
                    205:
                    206:   s = atan(ToReal(ARG0(arg)));
                    207:   MKReal(s,*rp);
                    208: }
                    209:
                    210: void Pdlog(NODE arg,Real *rp)
                    211: {
                    212:   double s;
                    213:
                    214:   s = log(ToReal(ARG0(arg)));
                    215:   MKReal(s,*rp);
                    216: }
                    217:
                    218: void Pdexp(NODE arg,Real *rp)
                    219: {
                    220:   double s;
                    221:
                    222:   s = exp(ToReal(ARG0(arg)));
                    223:   MKReal(s,*rp);
                    224: }
                    225:
                    226: void Pdfloor(NODE arg,Z *rp)
                    227: {
                    228:   double d;
                    229:   mpz_t z;
                    230:
                    231:   if ( !ARG0(arg) ) {
                    232:     *rp = 0;
                    233:     return;
                    234:   }
                    235:   d = floor(ToReal(ARG0(arg)));
                    236:   mpz_init(z);
                    237:   mpz_set_d(z,d);
                    238:   MPZTOZ(z,*rp);
                    239: }
                    240:
                    241: void Pdceil(NODE arg,Z *rp)
                    242: {
                    243:   double d;
                    244:   mpz_t z;
                    245:
                    246:   if ( !ARG0(arg) ) {
                    247:     *rp = 0;
                    248:     return;
                    249:   }
                    250:   d = ceil(ToReal(ARG0(arg)));
                    251:   mpz_init(z);
                    252:   mpz_set_d(z,d);
                    253:   MPZTOZ(z,*rp);
                    254: }
                    255:
                    256: void Pdrint(NODE arg,Z *rp)
                    257: {
                    258:   double d;
                    259:   mpz_t z;
                    260:
                    261:   if ( !ARG0(arg) ) {
                    262:     *rp = 0;
                    263:     return;
                    264:   }
                    265:   d = ToReal(ARG0(arg));
                    266: #if defined(VISUAL) || defined(__MINGW32__)
                    267:   d = d>0 ? floor(d+0.5) : ceil(d-0.5);
                    268: #else
                    269:   d = rint(d);
                    270: #endif
                    271:   mpz_init(z);
                    272:   mpz_set_d(z,d);
                    273:   MPZTOZ(z,*rp);
                    274: }
                    275:
                    276: void Pdisnan(NODE arg,Z *rp)
                    277: {
                    278:   Real r;
                    279:   double d;
                    280: #if defined(VISUAL) || defined(__MINGW32__)
                    281:     int c;
                    282: #endif
                    283:
                    284:   r = (Real)ARG0(arg);
                    285:   if ( !r || !NUM(r) || !REAL(r) ) {
                    286:     *rp = 0;
                    287:     return;
                    288:   }
                    289:   d = ToReal(r);
                    290: #if defined(VISUAL) || defined(__MINGW32__)
                    291:   c = _fpclass(d);
                    292:   if ( c == _FPCLASS_SNAN || c == _FPCLASS_QNAN ) *rp = ONE;
1.2     ! noro      293:   else if ( c == _FPCLASS_PINF || c == _FPCLASS_NINF ) STOZ(2,*rp);
1.1       noro      294: #else
                    295:   if ( isnan(d) ) *rp = ONE;
1.2     ! noro      296:   else if ( isinf(d) ) STOZ(2,*rp);
1.1       noro      297: #endif
                    298:   else *rp = 0;
                    299: }

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>